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