]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/GammaConv/AliAnalysisTaskPi0v2.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
Merge branch 'master' of https://git.cern.ch/reps/AliRoot
[u/mrichter/AliRoot.git] / PWGGA / GammaConv / AliAnalysisTaskPi0v2.cxx
CommitLineData
f7cf81bd 1#include <exception>
2#include "TRandom3.h"
3#include "TChain.h"
4#include "AliLog.h"
5#include "AliAnalysisTask.h"
6#include "AliAnalysisManager.h"
7#include "AliAnalysisTaskPi0v2.h"
8
9#include "AliESDEvent.h"
10#include "AliAODEvent.h"
11#include "AliCentrality.h"
12#include <iostream>
13
14#include "TFile.h"
15#include "AliOADBContainer.h"
16
e5b6e8a6 17#include "AliEPSelectionTask.h"
18
f7cf81bd 19// Author Daniel Lohner (Daniel.Lohner@cern.ch)
20
21using namespace std;
22
23ClassImp(AliAnalysisTaskPi0v2)
24
a91370e6 25const Double_t AliAnalysisTaskPi0v2::kGCfirstYBinSpectra = 0.;
26const Double_t AliAnalysisTaskPi0v2::kGClastYBinSpectra = 8.;
27
f7cf81bd 28//________________________________________________________________________
344100c4 29AliAnalysisTaskPi0v2::AliAnalysisTaskPi0v2(const char *name,Int_t harmonic) : AliAnalysisTaskSE(name),
30 fV0Reader(NULL),
31 fNCuts(0),
32 fConversionSelection(NULL),
33 fConversionGammas(NULL),
34 fNCentralityBins(knCentMax-1),
35 fCentrality(-1),
36 fCentralityBin(0),
37 fNBinsPhi(6),
38 fEP(NULL),
39 fUseTPCOnlyTracks(kTRUE),
40 fEtaMax(0.75),
41 fEtaGap(1),
42 fRPTPCEtaA(0),
43 fRPTPCEtaC(0),
44 fRPV0A(0),
45 fRPV0C(0),
46 fRPTPC(0),
47 fRPTPCEtaABF(0),
48 fRPTPCEtaCBF(0),
49 fRPV0ABF(0),
50 fRPV0CBF(0),
51 fRPTPCBF(0),
52 fEventCuts(NULL),
53 fConversionCuts(NULL),
54 fRandomizer(NULL),
55 fOutputList(NULL),
56 fMesonPDGCode(kPi0),
57 fDeltaPsiRP(0),
58 fRunNumber(0),
59 fRunIndex(0),
60 fNEPMethods(knEPMethod),
61 fFillQA(kTRUE),
62 fHarmonic(harmonic),
63 fPsiMax(2*TMath::Pi()/Double_t(harmonic)),
64 fPeriod("LHC10h"),
65 fIsAOD(kFALSE),
66 fSparseDist(NULL),
67 fHruns(NULL),
68 fDoEPFlattening(kTRUE),
69 fPeriodIndex(-1),
70 hNEvents(NULL),
71 hEventSelection(NULL),
72 hRPTPC(NULL),
73 hRPV0A(NULL),
74 hRPV0C(NULL),
75 hRPTPCAC(NULL),
76 hRPV0ATPC(NULL),
77 hRPV0CTPC(NULL),
78 hRPV0AC(NULL),
79 hCos2TPC(NULL),
80 hCos2V0ATPC(NULL),
81 hCos2V0CTPC(NULL),
82 hCos2V0AC(NULL),
83 hRPTPCEtaA(NULL),
84 hRPTPCEtaC(NULL),
85 hRPTPCEtaAC(NULL),
86 hCos2TPCEta(NULL),
87 hCos2V0ATPCEtaA(NULL),
88 hCos2V0ATPCEtaC(NULL),
89 hCos2V0CTPCEtaA(NULL),
90 hCos2V0CTPCEtaC(NULL),
91 hCos2SumWeights(NULL),
92 hEtaTPCEP(NULL),
93 hGammaMultCent(NULL),
94 hGammaPhi(NULL),
95 hMultChargedvsNGamma(NULL),
96 hMultChargedvsVZERO(NULL),
97 hMultChargedvsSPD(NULL),
98 hGammadNdPhi(NULL),
99 hGammaMultdPhiTRUE(NULL),
100 hGammaMultdPhiRECOTRUE(NULL),
101 hGammaMultTRUE(NULL),
102 hGammaMultRECOTRUE(NULL),
103 hGammaMultdPhi(NULL),
104 hGammaMult(NULL),
105 hGamma(NULL),
106 hGammaFull(NULL),
107 hCharged(NULL),
108 hPi0(NULL),
109 hPi0BG(NULL),
110
111 fMultV0(0x0),
112 fV0Cpol(0.),
113 fV0Apol(0.),
114 //hEPVertex(NULL)
115 hEPQA(NULL)
f7cf81bd 116{
344100c4 117 fInvMassRange[0]=0.05;
118 fInvMassRange[1]=0.3;
f7cf81bd 119
344100c4 120 for(Int_t ii=0;ii<knEPMethod;ii++)fEPSelectionMask[ii]=1;
e5b6e8a6 121
344100c4 122 fRandomizer=new TRandom3();
123 fRandomizer->SetSeed(0);
f7cf81bd 124
344100c4 125 for(Int_t i = 0; i < 4; ++i) {
126 fPhiDist[i] = 0;
127 }
e5b6e8a6 128
344100c4 129 for(Int_t ii=0;ii<knFlatPeriod;ii++){
130 for(Int_t jj=0;jj<knEP;jj++){
131 for(Int_t kk=0;kk<knCentMax;kk++){
132 fFlatc2[ii][jj][kk]=0;
133 fFlats2[ii][jj][kk]=0;
134 fFlatc4[ii][jj][kk]=0;
135 fFlats4[ii][jj][kk]=0;
136 }
137 }
33bf9167 138 }
33bf9167 139
344100c4 140 fCentralityBins[0]=0;
141 fCentralityBins[1]=5;
142 fCentralityBins[2]=10;
143 fCentralityBins[3]=20;
144 fCentralityBins[4]=30;
145 fCentralityBins[5]=40;
146 fCentralityBins[6]=50;
147 fCentralityBins[7]=60;
148 fCentralityBins[8]=70;
149 fCentralityBins[9]=80;
150
151 // Define input and output slots here
152 DefineInput(0, TChain::Class());
153 DefineOutput(1, TList::Class());
f7cf81bd 154}
155
156//________________________________________________________________________
157AliAnalysisTaskPi0v2::~AliAnalysisTaskPi0v2(){
158
344100c4 159 if(fRandomizer){
160 delete fRandomizer;
161 fRandomizer=NULL;
162 }
ca91a3e1 163
344100c4 164 if(fConversionSelection){
40a481fc 165 for(Int_t ii=0;ii<fNCuts;ii++)delete fConversionSelection[ii];
344100c4 166 delete[] fConversionSelection;
167 fConversionSelection=NULL;
168 }
70f90d01 169
344100c4 170 if(fEP){
171 delete fEP;
172 fEP=NULL;
e5b6e8a6 173 }
344100c4 174
175 if (fPeriod.CompareTo("LHC11h")==0){
176 for(Int_t i = 0; i < 4; i++) {
177 if(fPhiDist[i]){
178 delete fPhiDist[i];
179 fPhiDist[i] = NULL;
180 }
181 }
182 if(fHruns){
183 delete fHruns;
184 fHruns=NULL;
185 }
70f90d01 186 }
f7cf81bd 187}
188
189//________________________________________________________________________
190void AliAnalysisTaskPi0v2::UserCreateOutputObjects()
191{
344100c4 192 OpenFile(1);
f7cf81bd 193
344100c4 194 // GetConversionCuts
195 fConversionCuts=fV0Reader->GetConversionCuts();
196 fEventCuts=fV0Reader->GetEventCuts();
f7cf81bd 197
344100c4 198 // Flags
f7cf81bd 199
344100c4 200 Bool_t IsMC=AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler();
201 Bool_t IsHeavyIon=fEventCuts->IsHeavyIon();
f7cf81bd 202
344100c4 203 if(!IsHeavyIon||IsMC)fNEPMethods=1;
f7cf81bd 204
344100c4 205 // Create histograms
f7cf81bd 206
344100c4 207 if(fOutputList != NULL){
208 delete fOutputList;
209 fOutputList = NULL;
210 }
211 if(fOutputList == NULL){
212 fOutputList = new TList();
213 fOutputList->SetOwner(kTRUE);
214 }
f7cf81bd 215
344100c4 216 Int_t kGCnXBinsSpectra = Int_t((fInvMassRange[1]-fInvMassRange[0])*500); //500 for range 0 - 1
217 Double_t kGCfirstXBinSpectra = fInvMassRange[0];
218 Double_t kGClastXBinSpectra = fInvMassRange[1];
f7cf81bd 219
344100c4 220 Int_t nbinspi0[knbinsPi0]={kGCnYBinsSpectra,kGCnXBinsSpectra,fNBinsPhi,fNCentralityBins,fNEPMethods};
221 Double_t minpi0[knbinsPi0]={kGCfirstYBinSpectra,kGCfirstXBinSpectra,0,-0.5,-0.5};
222 Double_t maxpi0[knbinsPi0]={kGClastYBinSpectra,kGClastXBinSpectra,fPsiMax/2.,fNCentralityBins-0.5,fNEPMethods-0.5};
223 const char *binpi0[knbinsPi0]={"pt","mass","dPhi","centr","EPm"};
f7cf81bd 224
344100c4 225 Int_t nbinsg[knbinsGamma]={kGCnYBinsSpectra,fNBinsPhi,fNCentralityBins,fNEPMethods};
226 Double_t ming[knbinsGamma]={kGCfirstYBinSpectra,0,-0.5,-0.5};
227 Double_t maxg[knbinsGamma]={kGClastYBinSpectra,fPsiMax/2.,fNCentralityBins-0.5,fNEPMethods-0.5};
228 const char *bingamma[knbinsGamma]={"pt","dPhi","centr","EPm"};
f7cf81bd 229
344100c4 230 // Define Binning
f7cf81bd 231
344100c4 232 if(!IsMC){
f7cf81bd 233
344100c4 234 hPi0=new THnSparseF*[fNCuts+1];
235 hPi0BG=new THnSparseF*[fNCuts+1];
236 hGamma=new THnSparseF*[fNCuts+1];
f7cf81bd 237
344100c4 238 // Photon Cuts
239 for(Int_t iCut=0;iCut<fNCuts;iCut++){
f7cf81bd 240
344100c4 241 TList *fCutOutputList=new TList();
242 fCutOutputList->SetName(fConversionSelection[iCut]->GetCutString().Data());
243 fCutOutputList->SetOwner(kTRUE);
244 fOutputList->Add(fCutOutputList);
e5b6e8a6 245
344100c4 246 hPi0[iCut]=new THnSparseF("Pi0_Sparse","Pi0_Sparse",knbinsPi0,nbinspi0,minpi0,maxpi0);
247 for(Int_t i=0;i<knbinsPi0;i++) hPi0[iCut]->GetAxis(i)->SetName(binpi0[i]);
248 hPi0[iCut]->Sumw2();
249 fCutOutputList->Add(hPi0[iCut]);
e5b6e8a6 250
344100c4 251 hPi0BG[iCut]=new THnSparseF("Pi0BG_Sparse","Pi0BG_Sparse",knbinsPi0,nbinspi0,minpi0,maxpi0);
252 for(Int_t i=0;i<knbinsPi0;i++) hPi0BG[iCut]->GetAxis(i)->SetName(binpi0[i]);
253 hPi0BG[iCut]->Sumw2();
254 fCutOutputList->Add(hPi0BG[iCut]);
e5b6e8a6 255
344100c4 256 // Gamma
257 hGamma[iCut]=new THnSparseF("Gamma_Sparse","Gamma_Sparse",knbinsGamma,nbinsg,ming,maxg);
258 for(Int_t i=0;i<knbinsGamma;i++) hGamma[iCut]->GetAxis(i)->SetName(bingamma[i]);
259 hGamma[iCut]->Sumw2();
260 fCutOutputList->Add(hGamma[iCut]);
261 }
e5b6e8a6 262
344100c4 263 // no EP Flattening
264 Int_t iCut=0;
e5b6e8a6 265
344100c4 266 TList *fCutOutputList=new TList();
267 fCutOutputList->SetName(Form("%s_BF",fConversionSelection[iCut]->GetCutString().Data()));
268 fCutOutputList->SetOwner(kTRUE);
269 fOutputList->Add(fCutOutputList);
e5b6e8a6 270
344100c4 271 iCut=fNCuts;
e5b6e8a6 272
344100c4 273 hPi0[iCut]=new THnSparseF("Pi0_Sparse","Pi0_Sparse",knbinsPi0,nbinspi0,minpi0,maxpi0);
274 for(Int_t i=0;i<knbinsPi0;i++) hPi0[iCut]->GetAxis(i)->SetName(binpi0[i]);
275 hPi0[iCut]->Sumw2();
276 fCutOutputList->Add(hPi0[iCut]);
f7cf81bd 277
344100c4 278 hPi0BG[iCut]=new THnSparseF("Pi0BG_Sparse","Pi0BG_Sparse",knbinsPi0,nbinspi0,minpi0,maxpi0);
279 for(Int_t i=0;i<knbinsPi0;i++) hPi0BG[iCut]->GetAxis(i)->SetName(binpi0[i]);
280 hPi0BG[iCut]->Sumw2();
281 fCutOutputList->Add(hPi0BG[iCut]);
f7cf81bd 282
344100c4 283 // Gamma
284 hGamma[iCut]=new THnSparseF("Gamma_Sparse","Gamma_Sparse",knbinsGamma,nbinsg,ming,maxg);
285 for(Int_t i=0;i<knbinsGamma;i++) hGamma[iCut]->GetAxis(i)->SetName(bingamma[i]);
286 hGamma[iCut]->Sumw2();
287 fCutOutputList->Add(hGamma[iCut]);
f7cf81bd 288
f7cf81bd 289 }
290
344100c4 291 if(IsHeavyIon&&!IsMC){
292
293 // RP Calculation
294 TList *fRPList=new TList();
295 fRPList->SetName("Event Plane");
296 fRPList->SetOwner(kTRUE);
297 fOutputList->Add(fRPList);
298
299 hRPTPC=new TH2F("TPCAC" ,"TPC_AC" , fNCentralityBins,&fCentralityBins[0], 180, 0, fPsiMax);
300 hRPTPC->Sumw2();
301 fRPList->Add(hRPTPC);
302 hRPTPCEtaA=new TH2F("TPCEtaA" ,"TPC_EtaA" , fNCentralityBins,&fCentralityBins[0], 180, 0, fPsiMax);
303 hRPTPCEtaA->Sumw2();
304 fRPList->Add(hRPTPCEtaA);
305 hRPTPCEtaC=new TH2F("TPCEtaC" ,"TPC_EtaC" , fNCentralityBins,&fCentralityBins[0], 180, 0, fPsiMax);
306 hRPTPCEtaC->Sumw2();
307 fRPList->Add(hRPTPCEtaC);
308 hRPV0A=new TH2F("V0A" ,"VZERO_A" , fNCentralityBins,&fCentralityBins[0], 180, 0, fPsiMax);
309 hRPV0A->Sumw2();
310 fRPList->Add(hRPV0A);
311 hRPV0C=new TH2F("V0C" ,"VZERO_C" , fNCentralityBins,&fCentralityBins[0], 180, 0, fPsiMax);
312 hRPV0C->Sumw2();
313 fRPList->Add(hRPV0C);
314
315 hRPTPCAC=new TH2F("TPCA_TPCC" ,"TPCA_TPCC" , 180, 0, fPsiMax, 180, 0, fPsiMax);
316 hRPTPCAC->Sumw2();
317 fRPList->Add(hRPTPCAC);
318 hRPV0ATPC=new TH2F("V0A_TPC" ,"V0A_TPC" , 180, 0, fPsiMax, 180, 0, fPsiMax);
319 hRPV0ATPC->Sumw2();
320 fRPList->Add(hRPV0ATPC);
321 hRPV0CTPC=new TH2F("V0C_TPC" ,"V0C_TPC" , 180, 0, fPsiMax, 180, 0, fPsiMax);
322 hRPV0CTPC->Sumw2();
323 fRPList->Add(hRPV0CTPC);
324 hRPV0AC=new TH2F("V0A_V0C" ,"V0A_V0C" , 180, 0, fPsiMax, 180, 0, fPsiMax);
325 hRPV0AC->Sumw2();
326 fRPList->Add(hRPV0AC);
327 hRPTPCEtaAC=new TH2F("TPCEtaA_TPCEtaC" ,"TPCEtaA_TPCEtaC" , 180, 0, fPsiMax, 180, 0, fPsiMax);
328 hRPTPCEtaAC->Sumw2();
329 fRPList->Add(hRPTPCEtaAC);
330
331 Int_t nsystep=4;// 3 different weights + unflattened EP
332
333 hCos2TPC=new TH2F("Cos2_TPCAC" ,"Cos2_TPCAC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
334 hCos2TPC->Sumw2();
335 fRPList->Add(hCos2TPC);
336 hCos2TPCEta=new TH2F("Cos2_TPCEtaAC" ,"Cos2_TPCEtaAC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
337 hCos2TPCEta->Sumw2();
338 fRPList->Add(hCos2TPCEta);
339 hCos2V0ATPC=new TH2F("Cos2_V0ATPC" ,"Cos2_V0ATPC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
340 hCos2V0ATPC->Sumw2();
341 fRPList->Add(hCos2V0ATPC);
342 hCos2V0CTPC=new TH2F("Cos2_V0CTPC" ,"Cos2_V0CTPC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
343 hCos2V0CTPC->Sumw2();
344 fRPList->Add(hCos2V0CTPC);
345 hCos2V0AC=new TH2F("Cos2_V0AC" ,"Cos2_V0AC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
346 hCos2V0AC->Sumw2();
347 fRPList->Add(hCos2V0AC);
348 hCos2V0ATPCEtaA=new TH2F("Cos2_V0ATPCEtaA" ,"Cos2_V0ATPCEtaA" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
349 hCos2V0ATPCEtaA->Sumw2();
350 fRPList->Add(hCos2V0ATPCEtaA);
351 hCos2V0ATPCEtaC=new TH2F("Cos2_V0ATPCEtaC" ,"Cos2_V0ATPCEtaC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
352 hCos2V0ATPCEtaC->Sumw2();
353 fRPList->Add(hCos2V0ATPCEtaC);
354 hCos2V0CTPCEtaA=new TH2F("Cos2_V0CTPCEtaA" ,"Cos2_V0CTPCEtaA" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
355 hCos2V0CTPCEtaA->Sumw2();
356 fRPList->Add(hCos2V0CTPCEtaA);
357 hCos2V0CTPCEtaC=new TH2F("Cos2_V0CTPCEtaC" ,"Cos2_V0CTPCEtaC" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
358 hCos2V0CTPCEtaC->Sumw2();
359 fRPList->Add(hCos2V0CTPCEtaC);
360 hCos2SumWeights=new TH2F("Cos2_SumWeights" ,"Cos2_SumWeights" ,fNCentralityBins,&fCentralityBins[0],nsystep+1,-0.5,nsystep+0.5);
361 hCos2SumWeights->Sumw2();
362 fRPList->Add(hCos2SumWeights);
363
364 hEtaTPCEP=new TH2F("Eta_TPCEP" ,"Eta_TPCEP" ,fNCentralityBins,&fCentralityBins[0],100,-1,1);
365 hEtaTPCEP->Sumw2();
366 fRPList->Add(hEtaTPCEP);
367
368 /*const Int_t nepbins=4;
369 Int_t nbinsep[nepbins]={30,30,40,180};
370 Double_t minep[nepbins]={-0.015,0.17,-10,0};
371 Double_t maxep[nepbins]={0.015,0.20,10,fPsiMax};
372
373 hEPVertex=new THnSparseF("EP_Vertex","EP_Vertex",nepbins,nbinsep,minep,maxep);
374 fRPList->Add(hEPVertex);
375 */
376
377 const Int_t nRuns=275;
378 const Int_t nepbins=4;
379 Int_t nbinsep[nepbins]={fNCentralityBins,180,nRuns,5};
380 Double_t minep[nepbins]={-0.5,0,0.5,-0.5};
381 Double_t maxep[nepbins]={fNCentralityBins-0.5,fPsiMax,Double_t(nRuns)+0.5,4.5};
382 hEPQA=new THnSparseF("EP_QA","EP_QA",nepbins,nbinsep,minep,maxep);
383 fRPList->Add(hEPQA);
f7cf81bd 384 }
385
344100c4 386 TList *fPhotonQAList=new TList();
387 fPhotonQAList->SetName("Gamma_QA");
388 fPhotonQAList->SetOwner(kTRUE);
389 fOutputList->Add(fPhotonQAList);
390
391 if(fFillQA){
392 // Gamma QA
393 hGammaPhi=new TH2F*[fNCentralityBins];
394 for(Int_t m=0;m<fNCentralityBins;m++){
395 hGammaPhi[m]=new TH2F(Form("%d_GammaPhi",m),"GammaPhi",kGCnYBinsSpectra,kGCfirstYBinSpectra,kGClastYBinSpectra,360,0,2*TMath::Pi());
396 hGammaPhi[m]->Sumw2();
397 fPhotonQAList->Add(hGammaPhi[m]);
398 }
f7cf81bd 399
344100c4 400 hGammaMultCent=new TH2F("GammaMultvsCent","GammaMultvsCent",fNCentralityBins,&fCentralityBins[0], 60,-0.5,59.5);
401 hGammaMultCent->Sumw2();
402 fPhotonQAList->Add(hGammaMultCent);
403
404 hMultChargedvsSPD=new TH2F("Mult_ChargedvsSPD","Mult_ChargedvsSPD",250,0,2500, 250,0,5000);
405 hMultChargedvsSPD->Sumw2();
406 fPhotonQAList->Add(hMultChargedvsSPD);
407 hMultChargedvsVZERO=new TH2F("Mult_ChargedvsVZERO","Mult_ChargedvsVZERO",250,0,2500, 200,0,20000);
408 hMultChargedvsVZERO->Sumw2();
409 fPhotonQAList->Add(hMultChargedvsVZERO);
410 hMultChargedvsNGamma=new TH2F("Mult_ChargedvsNGamma","Mult_ChargedvsNGamma",250,0,2500,60,-0.5,59.5);
411 hMultChargedvsNGamma->Sumw2();
412 fPhotonQAList->Add(hMultChargedvsNGamma);
413
414 Int_t nbinsgmult[knbinsGammaMult]={kGCnYBinsSpectra,400,fNCentralityBins};
415 Double_t mingmult[knbinsGammaMult]={kGCfirstYBinSpectra,0,-0.5};
416 Double_t maxgmult[knbinsGammaMult]={kGClastYBinSpectra,8000,fNCentralityBins-0.5};
417 Double_t maxgmultdPhi[knbinsGammaMult]={kGClastYBinSpectra,2000,fNCentralityBins-0.5};
418 const char *bingammamult[knbinsGammaMult]={"pt","gmult","centr"};
419
420 if(AliAnalysisManager::GetAnalysisManager()->GetMCtruthEventHandler()){
421 hGammaMultdPhiTRUE=new THnSparseF("Gamma_MultdPhi_TRUE","Gamma_MultdPhi_TRUE",knbinsGammaMult,nbinsgmult,mingmult,maxgmultdPhi);
422 for(Int_t i=0;i<knbinsGammaMult;i++) hGammaMultdPhiTRUE->GetAxis(i)->SetName(bingammamult[i]);
423 hGammaMultdPhiTRUE->Sumw2();
424 fPhotonQAList->Add(hGammaMultdPhiTRUE);
425
426 hGammaMultTRUE=new THnSparseF("Gamma_Mult_TRUE","Gamma_Mult_TRUE",knbinsGammaMult,nbinsgmult,mingmult,maxgmult);
427 for(Int_t i=0;i<knbinsGammaMult;i++) hGammaMultTRUE->GetAxis(i)->SetName(bingammamult[i]);
428 hGammaMultTRUE->Sumw2();
429 fPhotonQAList->Add(hGammaMultTRUE);
430
431 hGammaMultdPhiRECOTRUE=new THnSparseF("Gamma_MultdPhi_RECOTRUE","Gamma_MultdPhi_RECOTRUE",knbinsGammaMult,nbinsgmult,mingmult,maxgmultdPhi);
432 for(Int_t i=0;i<knbinsGammaMult;i++) hGammaMultdPhiRECOTRUE->GetAxis(i)->SetName(bingammamult[i]);
433 hGammaMultdPhiRECOTRUE->Sumw2();
434 fPhotonQAList->Add(hGammaMultdPhiRECOTRUE);
435
436 hGammaMultRECOTRUE=new THnSparseF("Gamma_Mult_RECOTRUE","Gamma_Mult_RECOTRUE",knbinsGammaMult,nbinsgmult,mingmult,maxgmult);
437 for(Int_t i=0;i<knbinsGammaMult;i++) hGammaMultRECOTRUE->GetAxis(i)->SetName(bingammamult[i]);
438 hGammaMultRECOTRUE->Sumw2();
439 fPhotonQAList->Add(hGammaMultRECOTRUE);
440 }
f7cf81bd 441
344100c4 442 hGammaMultdPhi=new THnSparseF*[fNEPMethods];
443 hGammaMult=new THnSparseF*[fNEPMethods];
f7cf81bd 444
344100c4 445 hGammadNdPhi=new THnSparseF("Gamma_dNdPhi","Gamma_dNdPhi",knbinsGamma,nbinsg,ming,maxg);
446 for(Int_t i=0;i<knbinsGamma;i++) hGammadNdPhi->GetAxis(i)->SetName(bingamma[i]);
447 hGammadNdPhi->Sumw2();
448 fPhotonQAList->Add(hGammadNdPhi);
ca91a3e1 449
344100c4 450 for(Int_t iEP=0;iEP<fNEPMethods;iEP++){
451 hGammaMultdPhi[iEP]=new THnSparseF(Form("Gamma_MultdPhi_%d",iEP),"Gamma_MultdPhi",knbinsGammaMult,nbinsgmult,mingmult,maxgmultdPhi);
452 for(Int_t i=0;i<knbinsGammaMult;i++) hGammaMultdPhi[iEP]->GetAxis(i)->SetName(bingammamult[i]);
453 hGammaMultdPhi[iEP]->Sumw2();
454 fPhotonQAList->Add(hGammaMultdPhi[iEP]);
455
456 hGammaMult[iEP]=new THnSparseF(Form("Gamma_Mult_%d",iEP),"Gamma_Mult",knbinsGammaMult,nbinsgmult,mingmult,maxgmult);
457 for(Int_t i=0;i<knbinsGammaMult;i++) hGammaMult[iEP]->GetAxis(i)->SetName(bingammamult[i]);
458 hGammaMult[iEP]->Sumw2();
459 fPhotonQAList->Add(hGammaMult[iEP]);
460 }
461
462 const Int_t knbinsCharged=3;
463 Int_t nbinscharged[knbinsCharged]={fNBinsPhi,fNCentralityBins,fNEPMethods};
464 Double_t mincharged[knbinsCharged]={0,-0.5,-0.5};
465 Double_t maxcharged[knbinsCharged]={fPsiMax/2,fNCentralityBins-0.5,fNEPMethods-0.5};
466 hCharged=new THnSparseF("Charged","Charged",knbinsCharged,nbinscharged,mincharged,maxcharged);
467 hCharged->GetAxis(0)->SetName("dPhi");
468 hCharged->GetAxis(1)->SetName("centr");
469 hCharged->GetAxis(2)->SetName("EPm");
470 hCharged->Sumw2();
471 fPhotonQAList->Add(hCharged);
472
473 Int_t nbinsgfull[knbinsGamma]={kGCnYBinsSpectra,24,fNCentralityBins,fNEPMethods};
474 Double_t mingfull[knbinsGamma]={kGCfirstYBinSpectra,0,-0.5,-0.5};
475 Double_t maxgfull[knbinsGamma]={kGClastYBinSpectra,2*TMath::Pi(),fNCentralityBins-0.5,fNEPMethods-0.5};
476 hGammaFull=new THnSparseF("Gamma_Sparse_Full","Gamma_Sparse_Full",knbinsGamma,nbinsgfull,mingfull,maxgfull);
477 for(Int_t i=0;i<knbinsGamma;i++) hGammaFull->GetAxis(i)->SetName(bingamma[i]);
478 hGammaFull->Sumw2();
479 fPhotonQAList->Add(hGammaFull);
480 }
481
482 hNEvents=new TH1F("NEvents","NEvents",fNCentralityBins,&fCentralityBins[0]);
483 fPhotonQAList->Add(hNEvents);
484 hEventSelection=new TH1F("EventSelection","EventSelection",kEventSelected,0.5,kEventSelected+0.5);
485 hEventSelection->GetXaxis()->SetBinLabel(kEventIn,"in");
486 hEventSelection->GetXaxis()->SetBinLabel(kEventSelV0Reader,"v0reader");
487 hEventSelection->GetXaxis()->SetBinLabel(kEventCentrality,"centr");
488 hEventSelection->GetXaxis()->SetBinLabel(kEventRun,"run");
489 hEventSelection->GetXaxis()->SetBinLabel(kEventNoTPCEP,"no ep");
490 hEventSelection->GetXaxis()->SetBinLabel(kEventProcessEvent,"proc ev");
491 hEventSelection->GetXaxis()->SetBinLabel(kEventSelected,"selected");
492 fPhotonQAList->Add(hEventSelection);
493
494 // V0 Reader Cuts
495 TList *fV0ReaderCuts=fConversionCuts->GetCutHistograms();
496 fV0ReaderCuts->SetOwner(kTRUE);
497 fOutputList->Add(fV0ReaderCuts);
498
499 PostData(1, fOutputList);
f7cf81bd 500
501}
502
f7cf81bd 503//________________________________________________________________________
504Bool_t AliAnalysisTaskPi0v2::InitEvent(){
505
344100c4 506 if(!fV0Reader){AliError("Error: No V0 Reader and Pi0 Reconstructor");return kFALSE;}
507 if(!fV0Reader->IsEventSelected()){
508 hEventSelection->Fill(kEventSelV0Reader);
509 return kFALSE;
510 }
e5b6e8a6 511
344100c4 512 fConversionGammas=fV0Reader->GetReconstructedGammas();
f7cf81bd 513
344100c4 514 fIsAOD=(fInputEvent->IsA()==AliAODEvent::Class());
f7cf81bd 515
344100c4 516 if(!fConversionSelection){
517 AliError("No Cut Selection");
518 return kFALSE;
f7cf81bd 519 }
f7cf81bd 520
344100c4 521 if(!SetCentrality()){
522 hEventSelection->Fill(kEventCentrality);
523 return kFALSE;
a280ac15 524 }
a280ac15 525
344100c4 526 if(fEventCuts->IsHeavyIon()&&!fMCEvent){
f7cf81bd 527
344100c4 528 if(fRunNumber!=fInputEvent->GetRunNumber()){
529 fRunNumber=fInputEvent->GetRunNumber();
530 if (fRunNumber >= 136851 && fRunNumber <= 139515){fPeriod = "LHC10h";}
531 if (fRunNumber >= 166529 && fRunNumber <= 170593){fPeriod = "LHC11h";}
532 fRunIndex=GetRunIndex(fRunNumber);
533 fPeriodIndex=GetPeriodIndex(fPeriod);
534 LoadVZEROCalibration(fRunNumber); // Load Calibration for V0 Event Plane
535 LoadTPCCalibration(fRunNumber); // Load Calibration for TPC Event Plane
536 }
537 if(fRunIndex<0){
538 AliInfo("Run not selected");
539 hEventSelection->Fill(kEventRun);
540 return kFALSE;
541 }
e5b6e8a6 542
344100c4 543 // TPC Event Plane
544 if(!GetTPCEventPlane()){
545 hEventSelection->Fill(kEventNoTPCEP);
546 return kFALSE;
547 }
548 //fEP=fInputEvent->GetEventplane();
549 //if(!fEP)return kFALSE;
f7cf81bd 550
344100c4 551 fRPTPCBF=GetCorrectedTPCEPAngle(NULL,NULL,kFALSE);
552 fRPTPC=GetEventPlaneAngle(kTPC);
e5b6e8a6 553
344100c4 554 // TPC Eta Sub Events
555 fRPTPCEtaABF=GetTPCSubEPEta(kEPTPCEtaA);
556 fRPTPCEtaA=ApplyFlattening(fRPTPCEtaABF,kEPTPCEtaA);
557
558 fRPTPCEtaCBF=GetTPCSubEPEta(kEPTPCEtaC);
559 fRPTPCEtaC=ApplyFlattening(fRPTPCEtaCBF,kEPTPCEtaC);
560
561 // GetV0 Event Plane
562 GetV0EP(fInputEvent,fRPV0ABF,fRPV0CBF);
563 fRPV0A=ApplyFlattening(fRPV0ABF,kEPV0A);
564 fRPV0C=ApplyFlattening(fRPV0CBF,kEPV0C);
565 }
566 return kTRUE;
f7cf81bd 567}
568
e5b6e8a6 569
f7cf81bd 570//________________________________________________________________________
571void AliAnalysisTaskPi0v2::UserExec(Option_t *)
572{
344100c4 573 hEventSelection->Fill(kEventIn);
f7cf81bd 574
344100c4 575 if(!InitEvent())return;
f7cf81bd 576
344100c4 577 // Process Cuts
578 for(Int_t iCut=0;iCut<fNCuts;iCut++){
579 if(fConversionSelection[iCut]->ProcessEvent(fConversionGammas,fInputEvent,fMCEvent)){ // May only be called once per event!
580 if(!fMCEvent){
581 // Process EP methods
582 for(Int_t iEP=0;iEP<fNEPMethods;iEP++){
583 if(!fEPSelectionMask[iEP])continue; // dont fill THnSparse if not needed-> Save Memory
584 ProcessPi0s(iCut,EEventPlaneMethod(iEP));
585 ProcessGammas(iCut,EEventPlaneMethod(iEP));
586 }
587 }
588 // QA
589 if(fFillQA&&iCut==0)ProcessQA();
590 } else {
591 hEventSelection->Fill(kEventProcessEvent);
f7cf81bd 592 }
a280ac15 593 }
344100c4 594
595 // Fill N Events
596 hEventSelection->Fill(kEventSelected);
597 hNEvents->Fill(fCentrality);
598
599 // EventPlaneResolution
600 ProcessEventPlane();
601
602 PostData(1, fOutputList);
f7cf81bd 603}
604
605//________________________________________________________________________
606void AliAnalysisTaskPi0v2::ProcessPi0s(Int_t iCut,EEventPlaneMethod iEP){
607
344100c4 608 if(!fConversionSelection[iCut])return;
609 if(fConversionSelection[iCut]->GetNumberOfPhotons()==0)return;
f7cf81bd 610
344100c4 611 // Process Pi0s
f7cf81bd 612
344100c4 613 for(Int_t ii=0;ii<fConversionSelection[iCut]->GetNumberOfPi0s();ii++){
f7cf81bd 614
344100c4 615 AliAODConversionMother *pi0cand=fConversionSelection[iCut]->GetPi0(ii);
616 if(!pi0cand)continue;
f7cf81bd 617
344100c4 618 Double_t val[knbinsPi0];
619 val[kPi0Pt]=pi0cand->Pt();
620 val[kPi0Mass]=pi0cand->M();
621 val[kPi0dPhi]=GetPi0PhiwrtRP(pi0cand,iEP);
622 val[kPi0Cent]=fCentralityBin;
623 val[kPi0EPM]=Int_t(iEP);
624 hPi0[iCut]->Fill(val);
e5b6e8a6 625
344100c4 626 if(iCut==0){
627 // no flattening
628 val[kPi0dPhi]=GetPi0PhiwrtRP(pi0cand,iEP,kFALSE);
629 hPi0[fNCuts]->Fill(val);
630 }
e5b6e8a6 631 }
344100c4 632
633 // Pi0 BG
634 for(Int_t ii=0;ii<fConversionSelection[iCut]->GetNumberOfBGs();ii++){
635
636 AliAODConversionMother *pi0cand=fConversionSelection[iCut]->GetBG(ii);
637 if(!pi0cand)continue;
638
639 Double_t val[knbinsPi0];
640 val[kPi0Pt]=pi0cand->Pt();
641 val[kPi0Mass]=pi0cand->M();
642 val[kPi0dPhi]=GetPi0PhiwrtRP(pi0cand,iEP);
643 val[kPi0Cent]=fCentralityBin;
644 val[kPi0EPM]=Int_t(iEP);
645
646 hPi0BG[iCut]->Fill(val,pi0cand->GetWeight());
647
648 if(iCut==0){
649 // no flattening
650 val[kPi0dPhi]=GetPi0PhiwrtRP(pi0cand,iEP,kFALSE);
651 hPi0BG[fNCuts]->Fill(val);
652 }
e5b6e8a6 653 }
f7cf81bd 654}
655
656//________________________________________________________________________
657void AliAnalysisTaskPi0v2::ProcessGammas(Int_t iCut,EEventPlaneMethod iEP){
658
344100c4 659 if(!fConversionSelection[iCut]){
660 AliWarning("Conversion Selection does not exist");
661 return;
662 }
f7cf81bd 663
344100c4 664 for(Int_t ii=0;ii<fConversionSelection[iCut]->GetNumberOfPhotons();ii++){
e5b6e8a6 665
344100c4 666 AliAODConversionPhoton *gamma=fConversionSelection[iCut]->GetPhoton(ii);
667 Double_t val[knbinsGamma];
668 val[kGammaPt]=gamma->Pt();
669 val[kGammadPhi]=GetPhotonPhiwrtRP(gamma,iEP);
670 val[kGammaCent]=fCentralityBin;
671 val[kGammaEPM]=Int_t(iEP);
672 hGamma[iCut]->Fill(val);
f7cf81bd 673
344100c4 674 if(iCut==0){
675 // no flattening
676 val[kGammadPhi]=GetPhotonPhiwrtRP(gamma,iEP,kFALSE);
677 hGamma[fNCuts]->Fill(val);
678 }
ca91a3e1 679
344100c4 680 if(iCut==0&&fFillQA){
681 hGammadNdPhi->Fill(val);
682 Double_t EPAngle=GetEventPlaneAngle(iEP,gamma->Eta(),gamma,NULL);
683 Double_t dPhi=gamma->Phi()-EPAngle;
684 if(dPhi>=(2*TMath::Pi()))dPhi-=2*TMath::Pi();
685 if(dPhi<0)dPhi+=2*TMath::Pi();
686 val[kGammadPhi]=dPhi;
687 hGammaFull->Fill(val);
688 }
ca91a3e1 689 }
f7cf81bd 690}
691
692//________________________________________________________________________
693void AliAnalysisTaskPi0v2::ProcessQA(){
694
344100c4 695 if(!fConversionSelection[0]){
696 AliWarning("Conversion Selection does not exist");
697 return;
698 }
f7cf81bd 699
f7cf81bd 700
344100c4 701 AliStack *fMCStack=NULL;
702 if(fMCEvent)fMCStack=fMCEvent->Stack();
f7cf81bd 703
344100c4 704 // Multiplicity
705 Double_t multcharged=fConversionSelection[0]->GetNumberOfChargedTracks(fInputEvent);
706 Double_t multVZERO=fConversionSelection[0]->GetVZEROMult(fInputEvent);
707 Double_t multSPD=fConversionSelection[0]->GetSPDMult(fInputEvent);
f7cf81bd 708
344100c4 709 hMultChargedvsNGamma->Fill(multcharged,fConversionSelection[0]->GetNumberOfPhotons());
710 hMultChargedvsVZERO->Fill(multcharged,multVZERO);
711 hMultChargedvsSPD->Fill(multcharged,multSPD);
f7cf81bd 712
344100c4 713 // Efficiency Purity
714 Double_t valdPhi[knbinsGammaMult];
715 Double_t val[knbinsGammaMult];
f7cf81bd 716
344100c4 717 Int_t dNdPhi[fNBinsPhi];
718 Int_t ncharged;
f7cf81bd 719
344100c4 720 for(Int_t iEP=0;iEP<fNEPMethods;iEP++){
721 GetChargeddNdPhi(&dNdPhi[0],ncharged,iEP);
f7cf81bd 722
344100c4 723 // Reco
724 for(Int_t ii=0;ii<fConversionSelection[0]->GetNumberOfPhotons();ii++){
725 AliAODConversionPhoton *gamma=fConversionSelection[0]->GetPhoton(ii);
726 val[0]=gamma->Pt();
ca91a3e1 727 val[1]=ncharged;
728 val[2]=fCentralityBin;
344100c4 729
730 valdPhi[0]=gamma->Pt();
731 valdPhi[1]=dNdPhi[GetPhotonPhiBin(gamma,iEP)];
ca91a3e1 732 valdPhi[2]=fCentralityBin;
344100c4 733
734 hGammaMult[iEP]->Fill(val);
735 hGammaMultdPhi[iEP]->Fill(valdPhi);
736
737 // Gamma Phi
738 hGammaPhi[fCentralityBin]->Fill(gamma->Pt(),gamma->Phi());
739 hGammaMultCent->Fill(fCentrality,Float_t(fConversionSelection[0]->GetNumberOfPhotons()));
740
741 if(fMCStack){
742 if(gamma->IsTruePhoton(fMCStack)){
743 hGammaMultRECOTRUE->Fill(val);
744 hGammaMultdPhiRECOTRUE->Fill(valdPhi);
745 }
746 }
747 }
748
749 // MC Truth
750 if(fMCEvent){
751 for(Int_t i = 0; i < fMCStack->GetNprimary(); i++) {
752 TParticle* particle = (TParticle *)fMCStack->Particle(i);
753 if (!particle) continue;
754 if(fConversionCuts->PhotonIsSelectedMC(particle,fMCStack)){
755 TParticle *daughter=(TParticle *)fMCStack->Particle(particle->GetDaughter(0));
756 if(daughter){
757 val[0]=particle->Pt();
758 val[1]=ncharged;
759 val[2]=fCentralityBin;
760
761 valdPhi[0]=particle->Pt();
762 valdPhi[1]=dNdPhi[GetPhiBin(GetMCPhotonPhiwrtRP(particle,EEventPlaneMethod(iEP)))];
763 valdPhi[2]=fCentralityBin;
764
765 hGammaMultTRUE->Fill(val);
766 hGammaMultdPhiTRUE->Fill(valdPhi);
767 }
768 }
769 }
f7cf81bd 770 }
f7cf81bd 771 }
f7cf81bd 772}
773
774//________________________________________________________________________
775void AliAnalysisTaskPi0v2::GetPhotondNdPhi(Int_t *dNdPhi,Int_t iEP,Int_t iCut){
776
344100c4 777 for(Int_t iPhi=0;iPhi<fNBinsPhi;iPhi++)dNdPhi[iPhi]=0;
f7cf81bd 778
344100c4 779 for(Int_t ii=0;ii<fConversionSelection[iCut]->GetNumberOfPhotons();ii++){
780 AliAODConversionPhoton *gamma=fConversionSelection[iCut]->GetPhoton(ii);
781 Int_t phibin=GetPhotonPhiBin(gamma,iEP);
782 dNdPhi[phibin]++;
783 }
f7cf81bd 784}
785
786//________________________________________________________________________
787void AliAnalysisTaskPi0v2::GetChargeddNdPhi(Int_t *dNdPhi,Int_t &ntot,Int_t iEP){
788
344100c4 789 for(Int_t iPhi=0;iPhi<fNBinsPhi;iPhi++)dNdPhi[iPhi]=0;
790 ntot=0;
791
792 Double_t val[3];
793 val[1]=fCentralityBin;
794 val[2]=Int_t(iEP);
795
796 for(Int_t iTracks = 0; iTracks < fInputEvent->GetNumberOfTracks(); iTracks++){
797 AliVTrack* currentTrack = (AliVTrack*)fInputEvent->GetTrack(iTracks);
798 if(!currentTrack) continue;
799 if(TMath::Abs(currentTrack->Eta())>fEtaMax)continue;
800
801 Double_t phiwrt=GetChargedPhiwrtRP(currentTrack,EEventPlaneMethod(iEP));
802 Int_t phibin=GetPhiBin(phiwrt);
803
804 val[0]=phiwrt;
805 hCharged->Fill(val);
806
807 dNdPhi[phibin]++;
808 ntot++;
809 }
f7cf81bd 810}
811
812//________________________________________________________________________
813Int_t AliAnalysisTaskPi0v2::GetPhiBin(Double_t phiwrt){
344100c4 814 Double_t binrange=TMath::Pi()/(Double_t(fHarmonic*fNBinsPhi));
815 for(Int_t iPhi=0;iPhi<fNBinsPhi;iPhi++){
816 if(phiwrt>=(binrange*iPhi)&&phiwrt<(binrange*(iPhi+1)))return iPhi;
817 }
818 return -1;
f7cf81bd 819}
820
821//________________________________________________________________________
822Int_t AliAnalysisTaskPi0v2::GetPhotonPhiBin(AliAODConversionPhoton *gamma,Int_t iEP){
344100c4 823 Double_t phiwrt=GetPhotonPhiwrtRP(gamma,EEventPlaneMethod(iEP));
824 return GetPhiBin(phiwrt);
f7cf81bd 825}
826
827//________________________________________________________________________
e5b6e8a6 828Double_t AliAnalysisTaskPi0v2::GetPi0PhiwrtRP(AliAODConversionMother *pi0,EEventPlaneMethod iEP,Bool_t bDoFlattening){
f7cf81bd 829
344100c4 830 AliAODConversionPhoton *gamma0=dynamic_cast<AliAODConversionPhoton*>(fConversionGammas->At(pi0->GetLabel1()));
831 AliAODConversionPhoton *gamma1=dynamic_cast<AliAODConversionPhoton*>(fConversionGammas->At(pi0->GetLabel2()));
f7cf81bd 832
344100c4 833 Double_t EPAngle=GetEventPlaneAngle(iEP,pi0->Eta(),gamma0,gamma1,bDoFlattening);
f7cf81bd 834
344100c4 835 return GetPhiwrtRP(pi0->Phi()-EPAngle);
f7cf81bd 836}
837
838//________________________________________________________________________
e5b6e8a6 839Double_t AliAnalysisTaskPi0v2::GetPhotonPhiwrtRP(AliAODConversionPhoton *gamma,EEventPlaneMethod iEP,Bool_t bDoFlattening){
f7cf81bd 840
344100c4 841 Double_t EPAngle=GetEventPlaneAngle(iEP,gamma->Eta(),gamma,NULL,bDoFlattening);
f7cf81bd 842
344100c4 843 return GetPhiwrtRP(gamma->Phi()-EPAngle);
f7cf81bd 844}
845
846//________________________________________________________________________
e5b6e8a6 847Double_t AliAnalysisTaskPi0v2::GetMCPhotonPhiwrtRP(TParticle *gamma,EEventPlaneMethod iEP,Bool_t bDoFlattening){
f7cf81bd 848
344100c4 849 Double_t EPAngle=GetEventPlaneAngle(iEP,gamma->Eta(),NULL,NULL,bDoFlattening);
f7cf81bd 850
344100c4 851 return GetPhiwrtRP(gamma->Phi()-EPAngle);
f7cf81bd 852}
853//________________________________________________________________________
e5b6e8a6 854Double_t AliAnalysisTaskPi0v2::GetChargedPhiwrtRP(AliVTrack *track,EEventPlaneMethod iEP,Bool_t bDoFlattening){
f7cf81bd 855
344100c4 856 Double_t EPAngle=GetEventPlaneAngle(iEP,track->Eta(),NULL,NULL,bDoFlattening);
f7cf81bd 857
344100c4 858 return GetPhiwrtRP(track->Phi()-EPAngle);
e5b6e8a6 859}
860//________________________________________________________________________
861Double_t AliAnalysisTaskPi0v2::GetPhiwrtRP(Double_t dPhi){
344100c4 862 Double_t newdPhi=TMath::Abs(dPhi); // Cos is symmetric
863 while(newdPhi>fPsiMax/2)newdPhi=TMath::Abs(newdPhi-fPsiMax);
864 return newdPhi;
f7cf81bd 865}
866
867//________________________________________________________________________
868void AliAnalysisTaskPi0v2::Terminate(Option_t *)
869{
870
871}
872
873//________________________________________________________________________
874void AliAnalysisTaskPi0v2::ProcessEventPlane()
875{
344100c4 876 if(!fMCEvent&&fEventCuts->IsHeavyIon()){
877
878 /* Double_t val[4];
879 val[0]=fInputEvent->GetPrimaryVertex()->GetX();
880 val[1]=fInputEvent->GetPrimaryVertex()->GetY();
881 val[2]=fInputEvent->GetPrimaryVertex()->GetZ();
882 val[3]=GetEventPlaneAngle(kTPC);
883 hEPVertex->Fill(val); */
884
885 // Run by run monitoring (before flattening)
886
887 Double_t val[4];
888 val[0]=fCentralityBin;
889 val[2]=fRunIndex;
890
891 val[1]=fRPTPCBF;
892 val[3]=kEPTPC;
893 hEPQA->Fill(val);
894
895 val[1]=fRPTPCEtaABF;
896 val[3]=kEPTPCEtaA;
897 hEPQA->Fill(val);
898 val[1]=fRPTPCEtaCBF;
899 val[3]=kEPTPCEtaC;
900 hEPQA->Fill(val);
901
902 val[1]=fRPV0ABF;
903 val[3]=kEPV0A;
904 hEPQA->Fill(val);
905 val[1]=fRPV0CBF;
906 val[3]=kEPV0C;
907 hEPQA->Fill(val);
908
909 // After Flattening
910
911 // TPC EP
912 Double_t PsiRP1BF=fEP->GetQsub1()->Phi()/Double_t(fHarmonic);
913 Double_t PsiRP2BF=fEP->GetQsub2()->Phi()/Double_t(fHarmonic);
914 Double_t PsiRP1=ApplyFlattening(PsiRP1BF,kEPTPC);
915 Double_t PsiRP2=ApplyFlattening(PsiRP2BF,kEPTPC);
916
917 hRPTPC->Fill(fCentrality,fRPTPC);
918 hRPTPCAC->Fill(PsiRP1,PsiRP2);
919
920 // TPC Eta Gap
921 hRPTPCEtaA->Fill(fCentrality,fRPTPCEtaA);
922 hRPTPCEtaC->Fill(fCentrality,fRPTPCEtaC);
923 hRPTPCEtaAC->Fill(fRPTPCEtaA,fRPTPCEtaC);
924
925 // V0
926
927 hRPV0A->Fill(fCentrality,fRPV0A);
928 hRPV0C->Fill(fCentrality,fRPV0C);
929
930 hRPV0ATPC->Fill(fRPV0A,fRPTPC);
931 hRPV0CTPC->Fill(fRPV0C,fRPTPC);
932 hRPV0AC->Fill(fRPV0A,fRPV0C);
933
934 Double_t cos2V0ATPC=TMath::Cos(Double_t(fHarmonic)*(fRPTPC-fRPV0A));
935 Double_t cos2V0CTPC=TMath::Cos(Double_t(fHarmonic)*(fRPTPC-fRPV0C));
936 Double_t cos2V0AV0C=TMath::Cos(Double_t(fHarmonic)*(fRPV0C-fRPV0A));
937 Double_t cos2TPCEta=TMath::Cos(Double_t(fHarmonic)*(fRPTPCEtaA-fRPTPCEtaC));
938 Double_t cos2TPC=TMath::Cos(Double_t(fHarmonic)*(PsiRP1-PsiRP2));
939 Double_t cos2V0ATPCEtaA=TMath::Cos(Double_t(fHarmonic)*(fRPV0A-fRPTPCEtaA));
940 Double_t cos2V0CTPCEtaA=TMath::Cos(Double_t(fHarmonic)*(fRPV0C-fRPTPCEtaA));
941 Double_t cos2V0ATPCEtaC=TMath::Cos(Double_t(fHarmonic)*(fRPV0A-fRPTPCEtaC));
942 Double_t cos2V0CTPCEtaC=TMath::Cos(Double_t(fHarmonic)*(fRPV0C-fRPTPCEtaC));
943
944 const Int_t nfill=4;
945 Double_t weight[nfill];
946 weight[0]=1.;// Fill unweighted
947 weight[1]=Float_t(fConversionSelection[0]->GetNumberOfPhotons());// Weight with Photon Mult
948 weight[2]=Float_t(fConversionSelection[0]->GetNumberOfChargedTracks(fInputEvent)); // Weight with charged Track Mult
949 weight[3]=Float_t(fConversionSelection[0]->GetVZEROMult(fInputEvent)); // Weight with V0 mult
950
951 for(Int_t i=0;i<nfill;i++){
952
953 hCos2V0ATPC->Fill(fCentrality,i,cos2V0ATPC*weight[i]);
954 hCos2V0CTPC->Fill(fCentrality,i,cos2V0CTPC*weight[i]);
955 hCos2V0AC->Fill(fCentrality,i,cos2V0AV0C*weight[i]);
956 hCos2TPCEta->Fill(fCentrality,i,cos2TPCEta*weight[i]);
957 hCos2TPC->Fill(fCentrality,i,cos2TPC*weight[i]);
958 hCos2V0ATPCEtaA->Fill(fCentrality,i,cos2V0ATPCEtaA*weight[i]);
959 hCos2V0ATPCEtaC->Fill(fCentrality,i,cos2V0ATPCEtaC*weight[i]);
960 hCos2V0CTPCEtaA->Fill(fCentrality,i,cos2V0CTPCEtaA*weight[i]);
961 hCos2V0CTPCEtaC->Fill(fCentrality,i,cos2V0CTPCEtaC*weight[i]);
962
963 hCos2SumWeights->Fill(fCentrality,i,weight[i]);
964 }
ca91a3e1 965
344100c4 966 // Fill Resolution before EP Flattening
967 Double_t cos2V0ATPCBF=TMath::Cos(Double_t(fHarmonic)*(fRPTPCBF-fRPV0ABF));
968 Double_t cos2V0CTPCBF=TMath::Cos(Double_t(fHarmonic)*(fRPTPCBF-fRPV0CBF));
969 Double_t cos2V0AV0CBF=TMath::Cos(Double_t(fHarmonic)*(fRPV0CBF-fRPV0ABF));
970 Double_t cos2TPCEtaBF=TMath::Cos(Double_t(fHarmonic)*(fRPTPCEtaABF-fRPTPCEtaCBF));
971 Double_t cos2TPCBF=TMath::Cos(Double_t(fHarmonic)*(PsiRP1BF-PsiRP2BF));
972 Double_t cos2V0ATPCEtaABF=TMath::Cos(Double_t(fHarmonic)*(fRPV0ABF-fRPTPCEtaABF));
973 Double_t cos2V0CTPCEtaABF=TMath::Cos(Double_t(fHarmonic)*(fRPV0CBF-fRPTPCEtaABF));
974 Double_t cos2V0ATPCEtaCBF=TMath::Cos(Double_t(fHarmonic)*(fRPV0ABF-fRPTPCEtaCBF));
975 Double_t cos2V0CTPCEtaCBF=TMath::Cos(Double_t(fHarmonic)*(fRPV0CBF-fRPTPCEtaCBF));
976
977 hCos2V0ATPC->Fill(fCentrality,nfill,cos2V0ATPCBF);
978 hCos2V0CTPC->Fill(fCentrality,nfill,cos2V0CTPCBF);
979 hCos2V0AC->Fill(fCentrality,nfill,cos2V0AV0CBF);
980 hCos2TPCEta->Fill(fCentrality,nfill,cos2TPCEtaBF);
981 hCos2TPC->Fill(fCentrality,nfill,cos2TPCBF);
982 hCos2V0ATPCEtaA->Fill(fCentrality,nfill,cos2V0ATPCEtaABF);
983 hCos2V0ATPCEtaC->Fill(fCentrality,nfill,cos2V0ATPCEtaCBF);
984 hCos2V0CTPCEtaA->Fill(fCentrality,nfill,cos2V0CTPCEtaABF);
985 hCos2V0CTPCEtaC->Fill(fCentrality,nfill,cos2V0CTPCEtaCBF);
986
987 hCos2SumWeights->Fill(fCentrality,nfill);
ca91a3e1 988
344100c4 989 }
f7cf81bd 990}
991
992//________________________________________________________________________
993TVector2 AliAnalysisTaskPi0v2::GetEPContribution(AliAODConversionPhoton *gamma){
344100c4 994 TVector2 q;
995 for(Int_t ii=0;ii<2;ii++){
996 AliVTrack *fCurrentTrack=AliConversionPhotonCuts::GetTrack(fInputEvent,gamma->GetTrackLabel(ii));
997 TVector2 qtrack=GetContributionEP(fCurrentTrack);
998 q+=qtrack;
999 }
1000 return q;
f7cf81bd 1001}
1002
1003//________________________________________________________________________
e5b6e8a6 1004Double_t AliAnalysisTaskPi0v2::GetCorrectedTPCEPAngle(AliAODConversionPhoton *gamma0,AliAODConversionPhoton *gamma1,Bool_t bDoFlattening){
344100c4 1005 // Correct Event Plane for Dilepton Tracks
1006 TVector2 q0(*fEP->GetQVector());
1007 if(gamma0)q0-=GetEPContribution(gamma0);
1008 if(gamma1)q0-=GetEPContribution(gamma1);
1009 Double_t EPangle=GetPsiInRange(q0.Phi()/Double_t(fHarmonic));
1010 if(bDoFlattening)EPangle=ApplyFlattening(EPangle,kEPTPC);
1011
1012 return EPangle;
f7cf81bd 1013}
1014
1015//________________________________________________________________________
e5b6e8a6 1016Double_t AliAnalysisTaskPi0v2::GetTPCSubEPEta(EEventPlane ep){
1017
344100c4 1018 Double_t etamin,etamax;
1019 switch(ep){
1020 case kEPTPCEtaA:
1021 etamin=fEtaGap/2;
1022 etamax=1;
1023 break;
1024 case kEPTPCEtaC:
1025 etamin=-1;
1026 etamax=-fEtaGap/2;
1027 break;
1028 default:
1029 return 0;
1030 }
e5b6e8a6 1031
344100c4 1032 TVector2 q;
1033 for(Int_t ii=0;ii<fInputEvent->GetNumberOfTracks();ii++){
1034 AliVTrack *fCurrentTrack=dynamic_cast<AliVTrack*>(fInputEvent->GetTrack(ii));
1035 if(!fCurrentTrack)continue;
1036 if(fCurrentTrack->Eta()>=etamin&&fCurrentTrack->Eta()<=etamax){
1037 TVector2 qtrack=GetContributionEP(fCurrentTrack);
1038 q+=qtrack;
1039 }
f7cf81bd 1040 }
e5b6e8a6 1041
344100c4 1042 Double_t phi=GetPsiInRange(q.Phi()/Double_t(fHarmonic));
e5b6e8a6 1043
344100c4 1044 return phi;
f7cf81bd 1045}
1046
1047//________________________________________________________________________
e5b6e8a6 1048TVector2 AliAnalysisTaskPi0v2::GetContributionEP(AliVTrack *track){
1049
344100c4 1050 TVector2 q(0,0);
e5b6e8a6 1051
344100c4 1052 TArrayF *fQContributionX=fEP->GetQContributionXArray();
1053 TArrayF *fQContributionY=fEP->GetQContributionYArray();
e5b6e8a6 1054
344100c4 1055 Int_t trackID=track->GetID();
e5b6e8a6 1056
344100c4 1057 if(fIsAOD){
1058 if((trackID>-1)&&fUseTPCOnlyTracks)return q;
1059 if((trackID<0)&&!fUseTPCOnlyTracks)return q;
1060 if (fUseTPCOnlyTracks) trackID = trackID*(-1) - 1;
1061 }
e5b6e8a6 1062
344100c4 1063 q.Set(fQContributionX->GetAt(trackID),fQContributionY->GetAt(trackID));
e5b6e8a6 1064
344100c4 1065 return q;
e5b6e8a6 1066}
1067
1068//________________________________________________________________________
1069Double_t AliAnalysisTaskPi0v2::GetPsiInRange(Double_t phi){
1070
344100c4 1071 Double_t newphi=phi;
1072 while(newphi<0)newphi+=fPsiMax;
1073 while(newphi>fPsiMax)newphi-=fPsiMax;
1074 return newphi;
e5b6e8a6 1075}
1076
1077//________________________________________________________________________
1078Double_t AliAnalysisTaskPi0v2::GetEventPlaneAngle(EEventPlaneMethod EPmethod,Double_t eta,AliAODConversionPhoton *gamma0,AliAODConversionPhoton *gamma1,Bool_t bDoFlattening)
f7cf81bd 1079{
344100c4 1080 // If arguments are not null, the contribution of these photons is subtracted from the TPC EP
1081 if(fEventCuts->IsHeavyIon()){
1082 // For MC select random EP angle in order to avoid correlations due to azimuth dependent efficiencies (ITS holes)
1083 if(fMCEvent){
1084 return fRandomizer->Uniform(0,fPsiMax);
1085 }
1086 switch(EPmethod){
1087 case kTPC:
1088 return GetCorrectedTPCEPAngle(gamma0,gamma1,bDoFlattening);
1089 case kTPCEtaGap:
1090 // Use opposite EP
1091 if(bDoFlattening){
1092 if(eta<0)return fRPTPCEtaA;
1093 else return fRPTPCEtaC;
1094 } else{
1095 if(eta<0)return fRPTPCEtaABF;
1096 else return fRPTPCEtaCBF;
1097 }
1098 case kV0A:
1099 if(bDoFlattening)return fRPV0A;
1100 else return fRPV0ABF;
1101 case kV0C:
1102 if(bDoFlattening)return fRPV0C;
1103 else return fRPV0CBF;
1104 default:
1105 return 0;
1106 }
f7cf81bd 1107 }
1108
344100c4 1109 // NO EP in pp mode
1110 return 0;
f7cf81bd 1111}
1112
1113///________________________________________________________________________
1114Bool_t AliAnalysisTaskPi0v2::SetCentrality(){
1115
344100c4 1116 // Set centrality bin for current event
1117 if(!fEventCuts->IsHeavyIon()){
1118 fCentrality=0;
1119 fCentralityBin=0;
1120 return kTRUE;
1121 }
f7cf81bd 1122
344100c4 1123 fCentrality=fEventCuts->GetCentrality(fInputEvent);
1124 if(fNCentralityBins>1){
1125 for(fCentralityBin=0;fCentralityBin<fNCentralityBins;fCentralityBin++){
1126 if(fCentrality>=fCentralityBins[fCentralityBin]&&fCentrality<fCentralityBins[fCentralityBin+1])return kTRUE;
1127 }
1128 return kFALSE;
f7cf81bd 1129 }
344100c4 1130 fCentralityBin=0;
1131 return kTRUE;
f7cf81bd 1132}
1133
1134//________________________________________________________________________
1135void AliAnalysisTaskPi0v2::SetCentralityBins(Double_t *bins,Int_t nbins)
1136{
33bf9167 1137
344100c4 1138 if(nbins>=knCentMax){AliError("Not enough slots");return;}
f7cf81bd 1139
344100c4 1140 // Set Centrality bins for analysis
1141 fNCentralityBins=nbins;
1142 for(Int_t ii=0;ii<=fNCentralityBins;ii++){
1143 fCentralityBins[ii]=bins[ii];
1144 }
f7cf81bd 1145}
1146
1147//________________________________________________________________________
e5b6e8a6 1148void AliAnalysisTaskPi0v2::SetCuts(AliConversionSelection **conversionselection,Int_t ncuts){
f7cf81bd 1149
344100c4 1150 if(fConversionSelection){
1d9e6011 1151 for(Int_t ii=0;ii<fNCuts;ii++)delete fConversionSelection;
344100c4 1152 delete[] fConversionSelection;
1153 fConversionSelection=NULL;
1154 }
1155 fNCuts=ncuts;
1156 fConversionSelection=new AliConversionSelection*[fNCuts];
1157 for(Int_t ii=0;ii<fNCuts;ii++)fConversionSelection[ii]=new AliConversionSelection(*conversionselection[ii]);
f7cf81bd 1158}
1159
f7cf81bd 1160//___________________________________________________________________________
1161Int_t AliAnalysisTaskPi0v2::GetRunIndex(Int_t run){
1162
344100c4 1163 switch(run){
1164 //LHC11h (131 runs)
1165 case 167902 : return 140;
1166 case 167903 : return 141;
1167 case 167909 : return 142;
1168 case 167915 : return 143;
1169 case 167920 : return 144;
1170 case 167985 : return 145;
1171 case 167986 : return 146;
1172 case 167987 : return 147;
1173 case 167988 : return 148;
1174 case 168066 : return 149;
1175 case 168068 : return 150;
1176 case 168069 : return 151;
1177 case 168076 : return 152;
1178 case 168103 : return 153;
1179 case 168104 : return 154;
1180 case 168105 : return 155;
1181 case 168107 : return 156;
1182 case 168108 : return 157;
1183 case 168115 : return 158;
1184 case 168212 : return 159;
1185 case 168310 : return 160;
1186 case 168311 : return 161;
1187 case 168322 : return 162;
1188 case 168325 : return 163;
1189 case 168341 : return 164;
1190 case 168342 : return 165;
1191 case 168361 : return 166;
1192 case 168362 : return 167;
1193 case 168458 : return 168;
1194 case 168460 : return 169;
1195 case 168461 : return 170;
1196 case 168464 : return 171;
1197 case 168467 : return 172;
1198 case 168511 : return 173;
1199 case 168512 : return 174;
1200 case 168514 : return 175;
1201 case 168777 : return 176;
1202 case 168826 : return 177;
1203 case 168984 : return 178;
1204 case 168988 : return 179;
1205 case 168992 : return 180;
1206 case 169035 : return 181;
1207 case 169040 : return 182;
1208 case 169044 : return 183;
1209 case 169045 : return 184;
1210 case 169091 : return 185;
1211 case 169094 : return 186;
1212 case 169099 : return 187;
1213 case 169138 : return 188;
1214 case 169143 : return 189;
1215 case 169144 : return 190;
1216 case 169145 : return 191;
1217 case 169148 : return 192;
1218 case 169156 : return 193;
1219 case 169160 : return 194;
1220 case 169167 : return 195;
1221 case 169238 : return 196;
1222 case 169411 : return 197;
1223 case 169415 : return 198;
1224 case 169417 : return 199;
1225 case 169418 : return 200;
1226 case 169419 : return 201;
1227 case 169420 : return 202;
1228 case 169475 : return 203;
1229 case 169498 : return 204;
1230 case 169504 : return 205;
1231 case 169506 : return 206;
1232 case 169512 : return 207;
1233 case 169515 : return 208;
1234 case 169550 : return 209;
1235 case 169553 : return 210;
1236 case 169554 : return 211;
1237 case 169555 : return 212;
1238 case 169557 : return 213;
1239 case 169584 : return 214;
1240 case 169586 : return 215;
1241 case 169587 : return 216;
1242 case 169588 : return 217;
1243 case 169590 : return 218;
1244 case 169591 : return 219;
1245 case 169835 : return 220;
1246 case 169837 : return 221;
1247 case 169838 : return 222;
1248 case 169846 : return 223;
1249 case 169855 : return 224;
1250 case 169858 : return 225;
1251 case 169859 : return 226;
1252 case 169922 : return 227;
1253 case 169923 : return 228;
1254 case 169956 : return 229;
1255 case 169965 : return 230;
1256 case 169975 : return 231;
1257 case 169981 : return 232;
1258 case 170027 : return 233;
1259 case 170036 : return 234;
1260 case 170038 : return 235;
1261 case 170040 : return 236;
1262 case 170081 : return 237;
1263 case 170083 : return 238;
1264 case 170084 : return 239;
1265 case 170085 : return 240;
1266 case 170088 : return 241;
1267 case 170089 : return 242;
1268 case 170091 : return 243;
1269 case 170152 : return 244;
1270 case 170155 : return 245;
1271 case 170159 : return 246;
1272 case 170163 : return 247;
1273 case 170193 : return 248;
1274 case 170195 : return 249;
1275 case 170203 : return 250;
1276 case 170204 : return 251;
1277 case 170207 : return 252;
1278 case 170208 : return 253;
1279 case 170228 : return 254;
1280 case 170230 : return 255;
1281 case 170268 : return 256;
1282 case 170269 : return 257;
1283 case 170270 : return 258;
1284 case 170306 : return 259;
1285 case 170308 : return 260;
1286 case 170309 : return 261;
1287 case 170311 : return 262;
1288 case 170312 : return 263;
1289 case 170313 : return 264;
1290 case 170315 : return 265;
1291 case 170387 : return 266;
1292 case 170388 : return 267;
1293 case 170556 : return 268;
1294 case 170572 : return 269;
1295 case 170593 : return 270;
1296
1297 //LHC10h (137 runs)
1298 case 139517 : return 137;
1299 case 139514 : return 136;
1300 case 139513 : return 135;
1301 case 139511 : return 134;
1302 case 139510 : return 133;
1303 case 139507 : return 132;
1304 case 139505 : return 131;
1305 case 139504 : return 130;
1306 case 139503 : return 129;
1307 case 139470 : return 128;
1308 case 139467 : return 127;
1309 case 139466 : return 126;
1310 case 139465 : return 125;
1311 case 139440 : return 124;
1312 case 139439 : return 123;
1313 case 139438 : return 122;
1314 case 139437 : return 121;
1315 case 139360 : return 120;
1316 case 139329 : return 119;
1317 case 139328 : return 118;
1318 case 139314 : return 117;
1319 case 139311 : return 116;
1320 case 139310 : return 115;
1321 case 139309 : return 114;
1322 case 139308 : return 113;
1323 case 139173 : return 112;
1324 case 139172 : return 111;
1325 case 139110 : return 110;
1326 case 139107 : return 109;
1327 case 139105 : return 108;
1328 case 139104 : return 107;
1329 case 139042 : return 106;
1330 case 139038 : return 105;
1331 case 139037 : return 104;
1332 case 139036 : return 103;
1333 case 139029 : return 102;
1334 case 139028 : return 101;
1335 case 138983 : return 100;
1336 case 138982 : return 99;
1337 case 138980 : return 98;
1338 case 138979 : return 97;
1339 case 138978 : return 96;
1340 case 138977 : return 95;
1341 case 138976 : return 94;
1342 case 138973 : return 93;
1343 case 138972 : return 92;
1344 case 138965 : return 91;
1345 case 138924 : return 90;
1346 case 138872 : return 89;
1347 case 138871 : return 88;
1348 case 138870 : return 87;
1349 case 138837 : return 86;
1350 case 138830 : return 85;
1351 case 138828 : return 84;
1352 case 138826 : return 83;
1353 case 138796 : return 82;
1354 case 138795 : return 81;
1355 case 138742 : return 80;
1356 case 138732 : return 79;
1357 case 138730 : return 78;
1358 case 138666 : return 77;
1359 case 138662 : return 76;
1360 case 138653 : return 75;
1361 case 138652 : return 74;
1362 case 138638 : return 73;
1363 case 138624 : return 72;
1364 case 138621 : return 71;
1365 case 138583 : return 70;
1366 case 138582 : return 69;
1367 // case 138579 : return 68;
1368 // case 138578 : return 67;
1369 case 138534 : return 66;
1370 case 138469 : return 65;
1371 case 138442 : return 64;
1372 case 138439 : return 63;
1373 case 138438 : return 62;
1374 case 138396 : return 61;
1375 case 138364 : return 60;
1376 case 138359 : return 59;
1377 case 138275 : return 58;
1378 case 138225 : return 57;
1379 case 138201 : return 56;
1380 case 138200 : return 55;
1381 case 138197 : return 54;
1382 case 138192 : return 53;
1383 case 138190 : return 52;
1384 case 138154 : return 51;
1385 case 138153 : return 50;
1386 case 138151 : return 49;
1387 case 138150 : return 48;
1388 case 138126 : return 47;
1389 case 138125 : return 46;
1390 case 137848 : return 45;
1391 case 137847 : return 44;
1392 case 137844 : return 43;
1393 case 137843 : return 42;
1394 case 137752 : return 41;
1395 case 137751 : return 40;
1396 case 137748 : return 39;
1397 case 137724 : return 38;
1398 case 137722 : return 37;
1399 case 137718 : return 36;
1400 case 137704 : return 35;
1401 case 137693 : return 34;
1402 case 137692 : return 33;
1403 case 137691 : return 32;
1404 case 137689 : return 31;
1405 case 137686 : return 30;
1406 case 137685 : return 29;
1407 case 137639 : return 28;
1408 case 137638 : return 27;
1409 case 137608 : return 26;
1410 case 137595 : return 25;
1411 case 137549 : return 24;
1412 case 137546 : return 23;
1413 case 137544 : return 22;
1414 case 137541 : return 21;
1415 case 137539 : return 20;
1416 case 137531 : return 19;
1417 case 137530 : return 18;
1418 case 137443 : return 17;
1419 case 137441 : return 16;
1420 case 137440 : return 15;
1421 case 137439 : return 14;
1422 case 137434 : return 13;
1423 case 137432 : return 12;
1424 case 137431 : return 11;
1425 case 137430 : return 10;
1426 case 137366 : return 9;
1427 case 137243 : return 8;
1428 case 137236 : return 7;
1429 case 137235 : return 6;
1430 case 137232 : return 5;
1431 case 137231 : return 4;
1432 case 137165 : return 3;
1433 case 137162 : return 2;
1434 case 137161 : return 1;
1435
1436 // Default
1437 //default : return -1;
1438 default : return 275;
1439 }
f7cf81bd 1440}
1441
1442//____________________________________________________________________________
e5b6e8a6 1443void AliAnalysisTaskPi0v2::GetV0EP(AliVEvent * event,Double_t &rpv0a,Double_t &rpv0c){
f7cf81bd 1444
344100c4 1445 if (fPeriod.CompareTo("LHC10h")==0){
1446 // Corrected VZERO EP (from AliAnalysisTaskPi0Flow)
1447 //VZERO data
1448 AliESDVZERO* esdV0 = (AliESDVZERO*)event->GetVZEROData();
1449
1450 //reset Q vector info
1451 Double_t Qxa2 = 0, Qya2 = 0;
1452 Double_t Qxc2 = 0, Qyc2 = 0;
1453
1454 for (Int_t iv0 = 0; iv0 < 64; iv0++) {
1455 Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
1456 Float_t multv0 = esdV0->GetMultiplicity(iv0);
1457 Double_t lqx=TMath::Cos(fHarmonic*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
1458 Double_t lqy=TMath::Sin(fHarmonic*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
1459 if (iv0 < 32){ // V0C
1460 Qxc2 += lqx;
1461 Qyc2 += lqy;
1462 } else { // V0A
1463 Qxa2 += lqx;
1464 Qya2 += lqy;
1465 }
1466 }
f7cf81bd 1467
344100c4 1468 Int_t iC = -1;
1469 // centrality bins
1470 if(fCentrality < 5) iC = 0;
1471 else if(fCentrality < 10) iC = 1;
1472 else if(fCentrality < 20) iC = 2;
1473 else if(fCentrality < 30) iC = 3;
1474 else if(fCentrality < 40) iC = 4;
1475 else if(fCentrality < 50) iC = 5;
1476 else if(fCentrality < 60) iC = 6;
1477 else if(fCentrality < 70) iC = 7;
1478 else iC = 8;
1479
1480 //grab for each centrality the proper histo with the Qx and Qy to do the recentering
1481 Double_t Qxamean2 = fMeanQ[iC][1][0];
1482 Double_t Qxarms2 = fWidthQ[iC][1][0];
1483 Double_t Qyamean2 = fMeanQ[iC][1][1];
1484 Double_t Qyarms2 = fWidthQ[iC][1][1];
1485
1486 Double_t Qxcmean2 = fMeanQ[iC][0][0];
1487 Double_t Qxcrms2 = fWidthQ[iC][0][0];
1488 Double_t Qycmean2 = fMeanQ[iC][0][1];
1489 Double_t Qycrms2 = fWidthQ[iC][0][1];
1490
1491 Double_t QxaCor2 = (Qxa2 - Qxamean2)/Qxarms2;
1492 Double_t QyaCor2 = (Qya2 - Qyamean2)/Qyarms2;
1493 Double_t QxcCor2 = (Qxc2 - Qxcmean2)/Qxcrms2;
1494 Double_t QycCor2 = (Qyc2 - Qycmean2)/Qycrms2;
1495 rpv0a = TMath::ATan2(QyaCor2, QxaCor2)/Double_t(fHarmonic);
1496 rpv0c = TMath::ATan2(QycCor2, QxcCor2)/Double_t(fHarmonic);
1497
1498 //rpv0a = TMath::ATan2(Qya2, Qxa2)/Double_t(fHarmonic);
1499 //rpv0c = TMath::ATan2(Qyc2, Qxc2)/Double_t(fHarmonic);
1500
1501 // cout<<"Compare v"<<fHarmonic<<" "<<rpv0a<<" "<<fInputEvent->GetEventplane()->GetEventplane("V0A",fInputEvent,fHarmonic)<<endl;
1502 }
1503
1504 if (fPeriod.CompareTo("LHC11h")==0){
1505 AliEventplane *eventEP=fInputEvent->GetEventplane();
1506 Double_t qx,qy;
1507 rpv0a=eventEP->CalculateVZEROEventPlane(fInputEvent,8,fHarmonic,qx,qy);
1508 rpv0c=eventEP->CalculateVZEROEventPlane(fInputEvent,9,fHarmonic,qx,qy);
1509 }
f7cf81bd 1510
344100c4 1511 //AliEventplane *ep=fInputEvent->GetEventplane();
1512 //cout<<fHarmonic<<" A "<<ep->GetEventplane("V0A",fInputEvent,fHarmonic)<<" "<<rpv0a<<" C "<<ep->GetEventplane("V0C",fInputEvent,fHarmonic)<<" "<<rpv0c<<endl;
33bf9167 1513
344100c4 1514 rpv0a=GetPsiInRange(rpv0a);
1515 rpv0c=GetPsiInRange(rpv0c);
33bf9167 1516
344100c4 1517
f7cf81bd 1518}
1519
f7cf81bd 1520//_____________________________________________________________________________
1d9e6011 1521void AliAnalysisTaskPi0v2::LoadVZEROCalibration(Int_t run){
e5b6e8a6 1522
344100c4 1523 // VZERO Phi Weights and Recentering
1524 if (fPeriod.CompareTo("LHC10h")==0){
1525 TString oadbfilename = "$ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root";
1526 TFile *foadb = TFile::Open(oadbfilename.Data());
f7cf81bd 1527
344100c4 1528 if(!foadb){
1529 printf("OADB file %s cannot be opened\n",oadbfilename.Data());
1530 return;
1531 }
f7cf81bd 1532
344100c4 1533 AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr");
1534 if(!cont){
1535 printf("OADB object hMultV0BefCorr is not available in the file\n");
a280ac15 1536 return;
344100c4 1537 }
e5b6e8a6 1538
344100c4 1539 if(!(cont->GetObject(run))){
1540 printf("OADB object hMultV0BefCorr is not available for run %i (used run 137366)\n",run);
a280ac15 1541 run = 137366;
f7cf81bd 1542 }
344100c4 1543
1544 printf("Setting V0 calibration \n") ;
1545 fMultV0 = ((TH2F *) cont->GetObject(run))->ProfileX();
1546
1547 TF1 *fpol0 = new TF1("fpol0","pol0");
1548 fMultV0->Fit(fpol0,"0","",0,31);
1549 fV0Cpol = fpol0->GetParameter(0);
1550 fMultV0->Fit(fpol0,"0","",32,64);
1551 fV0Apol = fpol0->GetParameter(0);
1552
1553 for(Int_t iside=0;iside<2;iside++){
1554 for(Int_t icoord=0;icoord<2;icoord++){
1555 for(Int_t i=0;i < nCentrBinV0;i++){
1556 char namecont[100];
1557 if(iside==0 && icoord==0)
1558 snprintf(namecont,100,"hQxc%i_%i",fHarmonic,i);
1559 else if(iside==1 && icoord==0)
1560 snprintf(namecont,100,"hQxa%i_%i",fHarmonic,i);
1561 else if(iside==0 && icoord==1)
1562 snprintf(namecont,100,"hQyc%i_%i",fHarmonic,i);
1563 else if(iside==1 && icoord==1)
1564 snprintf(namecont,100,"hQya%i_%i",fHarmonic,i);
1565
1566 cont = (AliOADBContainer*) foadb->Get(namecont);
1567 if(!cont){
1568 printf("OADB object %s is not available in the file\n",namecont);
1569 return;
1570 }
1571
1572 if(!(cont->GetObject(run))){
1573 printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
1574 run = 137366;
1575 }
1576 fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
1577 fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
1578 }
1579 }
1580 }
e5b6e8a6 1581 }
1d9e6011 1582}
1583
1584//_____________________________________________________________________________
1585void AliAnalysisTaskPi0v2::LoadTPCCalibration(Int_t run){
f7cf81bd 1586
344100c4 1587 // TPC Event Plane Weights
1588 AliOADBContainer *fEPContainer=NULL;
1589 TString oadbfilename="";
f7cf81bd 1590
344100c4 1591 if (fPeriod.CompareTo("LHC10h")==0){
1592 // LHC10h
1593 if(fIsAOD){
1594 oadbfilename = (Form("%s/COMMON/EVENTPLANE/data/epphidist.aod.root", AliAnalysisManager::GetOADBPath()));
1595 } else{
1596 oadbfilename = (Form("%s/COMMON/EVENTPLANE/data/epphidist.root", AliAnalysisManager::GetOADBPath()));
1597 }
e5b6e8a6 1598
344100c4 1599 TFile foadb(oadbfilename);
1600 if(!foadb.IsOpen()) AliFatal(Form("Cannot open OADB file %s", oadbfilename.Data()));
1601
1602 AliInfo("Using Standard OADB");
1603 fEPContainer = (AliOADBContainer*) foadb.Get("epphidist");
1604 if (!fEPContainer) AliFatal("Cannot fetch OADB container for EP selection");
1605 foadb.Close();
1606
1607 fPhiDist[0] = (TH1F*) fEPContainer->GetObject(fRunNumber, "Default");
1608
1609 Bool_t emptybins;
1610
1611 int iter = 0;
1612 while (iter<3){
1613 emptybins = kFALSE;
1614 for (int i=1; i<fPhiDist[0]->GetNbinsX(); i++){
1615 if (!((fPhiDist[0]->GetBinContent(i))>0)) {
1616 emptybins = kTRUE;
1617 }
1618 }
1619 if (emptybins) {
1620 cout << "empty bins - rebinning!" << endl;
1621 fPhiDist[0]->Rebin();
1622 iter++;
1623 } else iter = 3;
f7cf81bd 1624 }
f7cf81bd 1625
344100c4 1626 if (emptybins) {
1627 AliError("After Maximum of rebinning still empty Phi-bins!!!");
1628 }
e5b6e8a6 1629 }
e5b6e8a6 1630
344100c4 1631 if (fPeriod.CompareTo("LHC11h")==0){
1632 // LHC11h
1633 oadbfilename = (Form("%s/COMMON/EVENTPLANE/data/epphidist2011.root", AliAnalysisManager::GetOADBPath()));
1634 TFile *foadb = TFile::Open(oadbfilename);
1635 if(!foadb->IsOpen()) AliFatal(Form("Cannot open OADB file %s", oadbfilename.Data()));
1636
1637 AliInfo("Using Standard OADB");
1638 fSparseDist = (THnSparse*) foadb->Get("Default");
1639 if (!fSparseDist) AliFatal("Cannot fetch OADB container for EP selection");
1640 foadb->Close();
1641 if(!fHruns){
1642 fHruns = (TH1F*)fSparseDist->Projection(0); //projection on run axis;
1643 fHruns->SetName("runsHisto");
1644 }
e5b6e8a6 1645
344100c4 1646 Int_t runbin=fHruns->FindBin(fRunNumber);
1647 if (fHruns->GetBinContent(runbin) > 1){
1648 fSparseDist->GetAxis(0)->SetRange(runbin,runbin);
1649 } else if(fHruns->GetBinContent(runbin) < 2){
1650 fSparseDist->GetAxis(0)->SetRange(1,2901); // not calibrated run, use integrated phi-weights
1651 AliInfo("Using integrated Phi-weights for this run");
1652 }
1653 for (Int_t i = 0; i<4 ;i++) {
1654 if(fPhiDist[i]){
1655 delete fPhiDist[i];
1656 fPhiDist[i] = 0x0;
1657 }
1658 if(i == 0){
1659 fSparseDist->GetAxis(1)->SetRange(1,1); // neg charge
1660 fSparseDist->GetAxis(2)->SetRange(1,1);} // neg eta
1661 if(i == 1){
1662 fSparseDist->GetAxis(1)->SetRange(2,2); // pos charge
1663 fSparseDist->GetAxis(2)->SetRange(1,1);} // neg eta
1664 if(i == 2){
1665 fSparseDist->GetAxis(1)->SetRange(1,1); // neg charge
1666 fSparseDist->GetAxis(2)->SetRange(2,2);} // pos eta
1667 if(i == 3){
1668 fSparseDist->GetAxis(1)->SetRange(2,2); // pos charge
1669 fSparseDist->GetAxis(2)->SetRange(2,2);} // pos eta
1670 fPhiDist[i] = (TH1F*)fSparseDist->Projection(3); // Projection on Phi
1671 fPhiDist[i]->SetName(Form("phidist%d%d",i,fRunNumber));
1672 fSparseDist->GetAxis(1)->SetRange(1,2); // reset axes
1673 fSparseDist->GetAxis(2)->SetRange(1,2);
1674 }
1675 fSparseDist->GetAxis(0)->SetRange(1,2901);// reset run axis
e5b6e8a6 1676 }
344100c4 1677 if (!fPhiDist[0]) AliFatal(Form("Cannot find OADB phi distribution for run %d", run));
e5b6e8a6 1678
1679}
1680
1681//_________________________________________________________________________
1682Int_t AliAnalysisTaskPi0v2::GetAODEPTrackFilterBit(){
1683
344100c4 1684 if(fUseTPCOnlyTracks){
1685 return 128;// TPC only with vertex constraint
1686 }
1687 return 1;// Use Global Tracks
e5b6e8a6 1688}
1689
1690//_________________________________________________________________________
1691TObjArray* AliAnalysisTaskPi0v2::GetEventPlaneTracks(Int_t &maxID)
1692{
344100c4 1693 TObjArray *tracklist=NULL;
1694 AliESDtrackCuts *fEPESDtrackCuts=AliESDtrackCuts::GetStandardTPCOnlyTrackCuts();
1695 fEPESDtrackCuts->SetPtRange(0.15,20.);
1696 fEPESDtrackCuts->SetEtaRange(-0.8,0.8);
1697 Int_t fAODfilterbit=GetAODEPTrackFilterBit();
1698
1699 if(fInputEvent->IsA()==AliESDEvent::Class()){
1700 maxID=fInputEvent->GetNumberOfTracks();
1701 tracklist=fEPESDtrackCuts->GetAcceptedTracks((AliESDEvent*)fInputEvent,fUseTPCOnlyTracks);
e5b6e8a6 1702 }
344100c4 1703
1704 if(fInputEvent->IsA()==AliAODEvent::Class()){
1705 // From AliEPSelectionTask
1706 tracklist = new TObjArray();
1707
1708 AliAODTrack *tr = 0;
1709 Int_t maxid1 = 0;
1710 Int_t maxidtemp = -1;
1711 Float_t ptlow = 0;
1712 Float_t ptup = 0;
1713 Float_t etalow = 0;
1714 Float_t etaup = 0;
1715 fEPESDtrackCuts->GetPtRange(ptlow,ptup);
1716 fEPESDtrackCuts->GetEtaRange(etalow,etaup);
1717 Int_t ntpc = fEPESDtrackCuts->GetMinNClusterTPC();
1718
1719 for (Int_t i = 0; i < fInputEvent->GetNumberOfTracks() ; i++){
1720 tr = (AliAODTrack*)fInputEvent->GetTrack(i);
1721 maxidtemp = tr->GetID();
1722 if(maxidtemp < 0 && fAODfilterbit != 128) continue;// id<0 means filter bit 128
1723 if(maxidtemp > -1 && fAODfilterbit == 128) continue;// id>01 means filter bit 1
1724 if (fAODfilterbit == 128) maxidtemp = maxidtemp*(-1) - 1;
1725 if (maxidtemp > maxid1) maxid1 = maxidtemp;
1726 if(tr->TestFilterBit(fAODfilterbit) && tr->Pt() < ptup && tr->Pt() > ptlow && tr->Eta() < etaup && tr->Eta() > etalow && tr->GetTPCNcls() > ntpc){
1727 tracklist->Add(tr);
1728 }
1729 }
1730 maxID = maxid1;
1731 }
1732 delete fEPESDtrackCuts;
1733 if(!tracklist)AliError("No tracklist");
1734 return tracklist;
f7cf81bd 1735}
1736
1737//____________________________________________________________________________
e5b6e8a6 1738Bool_t AliAnalysisTaskPi0v2::GetTPCEventPlane(){
1739
344100c4 1740 if(fEP){
1741 delete fEP;
1742 fEP=NULL;
1743 }
1744 fEP=new AliEventplane();
1745
1746 float mQx=0, mQy=0;
1747 float mQx1=0, mQy1=0, mQx2=0, mQy2=0;
1748 AliVTrack* track;
1749 Double_t weight;
1750 Int_t idtemp = -1;
1751 int trackcounter1=0, trackcounter2=0;
1752
1753 Int_t maxID=0;
1754
1755 TObjArray *tracklist=GetEventPlaneTracks(maxID);
1756
1757 fEP->GetQContributionXArray()->Set(maxID);
1758 fEP->GetQContributionYArray()->Set(maxID);
1759 fEP->GetQContributionXArraysub1()->Set(maxID);
1760 fEP->GetQContributionYArraysub1()->Set(maxID);
1761 fEP->GetQContributionXArraysub2()->Set(maxID);
1762 fEP->GetQContributionYArraysub2()->Set(maxID);
1763
1764 int nt = tracklist->GetEntries();
1765 for (int i=0; i<nt; i++){
1766 weight = 1;
1767 track = dynamic_cast<AliVTrack*> (tracklist->At(i));
1768 if (track) {
1769 // Fill Eta Distribution
1770 hEtaTPCEP->Fill(fCentrality,track->Eta());
1771
1772 weight=GetWeight(track);
1773 idtemp = track->GetID();
1774 // TPC only tracks have negative id ((-1)*IDESD - 1) in AOD
1775 if (fIsAOD && (fUseTPCOnlyTracks)) idtemp = idtemp*(-1) - 1;
1776
1777 Double_t qx=weight*cos(Double_t(fHarmonic)*track->Phi());
1778 Double_t qy=weight*sin(Double_t(fHarmonic)*track->Phi());
1779 fEP->GetQContributionXArray()->AddAt(qx,idtemp);
1780 fEP->GetQContributionYArray()->AddAt(qy,idtemp);
1781
1782 mQx += (qx);
1783 mQy += (qy);
1784
1785 // This loop splits the track set into 2 random subsets
1786 if( trackcounter1 < int(nt/2.) && trackcounter2 < int(nt/2.)){
1787 float random = fRandomizer->Rndm();
1788 if(random < .5){
1789 mQx1 += (qx);
1790 mQy1 += (qy);
1791 fEP->GetQContributionXArraysub1()->AddAt(qx,idtemp);
1792 fEP->GetQContributionYArraysub1()->AddAt(qy,idtemp);
1793 trackcounter1++;
1794 } else {
1795 mQx2 += (qx);
1796 mQy2 += (qy);
1797 fEP->GetQContributionXArraysub2()->AddAt(qx,idtemp);
1798 fEP->GetQContributionYArraysub2()->AddAt(qy,idtemp);
1799 trackcounter2++;
1800 }
1801 } else{
1802 if( trackcounter1 >= int(nt/2.)){
1803 mQx2 += (qx);
1804 mQy2 += (qy);
1805 fEP->GetQContributionXArraysub2()->AddAt(qx,idtemp);
1806 fEP->GetQContributionYArraysub2()->AddAt(qy,idtemp);
1807 trackcounter2++;
1808 } else {
1809 mQx1 += (qx);
1810 mQy1 += (qy);
1811 fEP->GetQContributionXArraysub1()->AddAt(qx,idtemp);
1812 fEP->GetQContributionYArraysub1()->AddAt(qy,idtemp);
1813 trackcounter1++;
1814 }
1815 }
e5b6e8a6 1816 }
e5b6e8a6 1817 }
f7cf81bd 1818
344100c4 1819 tracklist->Clear();
1820 delete tracklist;
1821 tracklist = NULL;
1822
1823 TVector2 *mQ=new TVector2();
1824 mQ->Set(mQx,mQy);
1825 Double_t EPAngle=mQ->Phi()/Double_t(fHarmonic);
f7cf81bd 1826
344100c4 1827 TVector2 *fQsub1=new TVector2();
1828 TVector2 *fQsub2=new TVector2();
1829 fQsub1->Set(mQx1,mQy1);
1830 fQsub2->Set(mQx2,mQy2);
f7cf81bd 1831
344100c4 1832 fEP->SetQVector(mQ);
1833 fEP->SetEventplaneQ(EPAngle);
1834 fEP->SetQsub(fQsub1,fQsub2);
1835 fEP->SetQsubRes(fQsub1->Phi()/Double_t(fHarmonic) - fQsub2->Phi()/Double_t(fHarmonic));
f7cf81bd 1836
344100c4 1837 Int_t ntracks=trackcounter1+trackcounter2;
1d9e6011 1838
344100c4 1839 //AliEventplane *ep=fInputEvent->GetEventplane();
1840 //if(fHarmonic==2)cout<<ep->GetQVector()->Phi()<<" "<<fEP->GetQVector()->Phi()<<endl;
33bf9167 1841
344100c4 1842 if(ntracks<3)return kFALSE;// <3 -> no subevents
1843 return kTRUE;
f7cf81bd 1844}
e5b6e8a6 1845
f7cf81bd 1846//____________________________________________________________________________
33bf9167 1847void AliAnalysisTaskPi0v2::SetFlatteningCoeff(EEventPlane ep,Int_t period,Int_t nCent,Double_t *cc2,Double_t *cs2,Double_t *cc4,Double_t *cs4){
e5b6e8a6 1848
344100c4 1849 if(nCent>knCentMax){AliError("Exceeds available Slots");return;}
1850 if(ep>=knEP||ep<0){AliError("Unknown EPMethod");return;}
1851 if(period>=knFlatPeriod||period<0){AliError("Unknown EPMethod");return;}
e5b6e8a6 1852
344100c4 1853 for(Int_t ic=0;ic<nCent;ic++){
1854 fFlatc2[period][ep][ic]=cc2[ic];
1855 fFlats2[period][ep][ic]=cs2[ic];
1856 fFlatc4[period][ep][ic]=cc4[ic];
1857 fFlats4[period][ep][ic]=cs4[ic];
1858 }
33bf9167 1859}
e5b6e8a6 1860
33bf9167 1861//____________________________________________________________________________
1862Int_t AliAnalysisTaskPi0v2::GetPeriodIndex(TString period){
e5b6e8a6 1863
344100c4 1864 if(period.CompareTo("LHC10h")==0)return 0;
1865 if(period.CompareTo("LHC11h")==0)return 1;
e5b6e8a6 1866
344100c4 1867 return -1;
33bf9167 1868}
e5b6e8a6 1869
33bf9167 1870//____________________________________________________________________________
1871Double_t AliAnalysisTaskPi0v2::ApplyFlattening(Double_t phi,EEventPlane ep){
e5b6e8a6 1872
344100c4 1873 Double_t newphi=phi;
e5b6e8a6 1874
344100c4 1875 if(fPeriodIndex>=0){
1876 Double_t c2=fFlatc2[fPeriodIndex][Int_t(ep)][fCentralityBin];
1877 Double_t s2=fFlats2[fPeriodIndex][Int_t(ep)][fCentralityBin];
1878 Double_t c4=fFlatc4[fPeriodIndex][Int_t(ep)][fCentralityBin];
1879 Double_t s4=fFlats4[fPeriodIndex][Int_t(ep)][fCentralityBin];
e5b6e8a6 1880
344100c4 1881 // Do correction
1882 newphi+=1/Double_t(fHarmonic)*(2*c2*sin(Double_t(fHarmonic)*phi)-2*s2*cos(Double_t(fHarmonic)*phi)+c4*sin(2*Double_t(fHarmonic)*phi)-s4*cos(2*Double_t(fHarmonic)*phi));
1883 newphi=GetPsiInRange(newphi);
e5b6e8a6 1884
344100c4 1885 }
e5b6e8a6 1886
344100c4 1887 return newphi;
e5b6e8a6 1888}
1889
1890//________________________________________________________________________
1891Double_t AliAnalysisTaskPi0v2::GetWeight(TObject* track1)
1892{
344100c4 1893 Double_t ptweight=1;
1894 AliVTrack* track = dynamic_cast<AliVTrack*>(track1);
1895 if (track) {
1896 if (track->Pt()<2) ptweight=track->Pt();
1897 else ptweight=2;
1898 }
1899 return ptweight*GetPhiWeight(track);
e5b6e8a6 1900}
f7cf81bd 1901
e5b6e8a6 1902//________________________________________________________________________
1903Double_t AliAnalysisTaskPi0v2::GetPhiWeight(TObject* track1)
1904{
344100c4 1905 Double_t phiweight=1;
1906 AliVTrack* track = dynamic_cast<AliVTrack*>(track1);
f7cf81bd 1907
344100c4 1908 TH1F *phiDist = 0x0;
1909 if(track) phiDist = SelectPhiDist(track);
e5b6e8a6 1910
344100c4 1911 if (phiDist && track) {
1912 Double_t nParticles = phiDist->Integral();
1913 Double_t nPhibins = phiDist->GetNbinsX();
e5b6e8a6 1914
344100c4 1915 Double_t Phi = track->Phi();
e5b6e8a6 1916
344100c4 1917 while (Phi<0) Phi += TMath::TwoPi();
1918 while (Phi>TMath::TwoPi()) Phi -= TMath::TwoPi();
e5b6e8a6 1919
344100c4 1920 Double_t PhiDistValue = phiDist->GetBinContent(1+TMath::FloorNint((track->Phi())*nPhibins/TMath::TwoPi()));
e5b6e8a6 1921
344100c4 1922 if (PhiDistValue > 0) phiweight = nParticles/nPhibins/PhiDistValue;
1923 }
1924 return phiweight;
f7cf81bd 1925}
e5b6e8a6 1926
1927//_________________________________________________________________________
1928TH1F* AliAnalysisTaskPi0v2::SelectPhiDist(AliVTrack *track)
1929{
344100c4 1930 if (fPeriod.CompareTo("LHC10h")==0) return fPhiDist[0];
1931 else if(fPeriod.CompareTo("LHC11h")==0){
1932 if (track->Charge() < 0){
1933 if(track->Eta() < 0.) return fPhiDist[0];
1934 else if (track->Eta() > 0.) return fPhiDist[2];
1935 }
1936 else if (track->Charge() > 0){
1937 if(track->Eta() < 0.) return fPhiDist[1];
1938 else if (track->Eta() > 0.) return fPhiDist[3];
1939 }
1940 }
1941 return 0;
e5b6e8a6 1942}