]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWG/FLOW/Tasks/AliAnalysisTaskVnV0.cxx
force compute MC label
[u/mrichter/AliRoot.git] / PWG / FLOW / Tasks / AliAnalysisTaskVnV0.cxx
CommitLineData
afa8df58 1#include "AliAnalysisTaskVnV0.h"
2
3// ROOT includes
4#include <TMath.h>
5
6// AliRoot includes
7#include "AliInputEventHandler.h"
8#include "AliAODEvent.h"
9#include "AliAODVertex.h"
10#include "AliAODTrack.h"
11#include "AliCentrality.h"
12#include "AliVHeader.h"
13#include "AliAODVZERO.h"
14#include "TFile.h"
15#include "AliOADBContainer.h"
16#include "TH2F.h"
17#include "TF1.h"
587d006a 18#include "AliGenHijingEventHeader.h"
19#include "AliMCEvent.h"
20#include "AliAODMCHeader.h"
21#include "AliAODMCParticle.h"
243fbce7 22#include "TChain.h"
afa8df58 23
24// STL includes
25//#include <iostream>
26//using namespace std;
27
28ClassImp(AliAnalysisTaskVnV0)
29
afa8df58 30//_____________________________________________________________________________
31AliAnalysisTaskVnV0::AliAnalysisTaskVnV0():
32 AliAnalysisTaskSE(),
587d006a 33 fVtxCut(10.0), // cut on |vertex| < fVtxCut
afa8df58 34 fEtaCut(0.8), // cut on |eta| < fEtaCut
35 fMinPt(0.15), // cut on pt > fMinPt
243fbce7 36 fV2(kTRUE),
37 fV3(kTRUE),
38 fIsMC(kFALSE),
39 fQAsw(kFALSE),
afa8df58 40 fRun(-1),
41 fList(new TList()),
42 fList2(new TList()),
587d006a 43 fList3(new TList()),
44 fList4(new TList()),
45 fMultV0(NULL),
afa8df58 46 fV0Cpol(100),
47 fV0Apol(100),
587d006a 48 fHResTPCv0A2(NULL),
49 fHResTPCv0C2(NULL),
50 fHResv0Cv0A2(NULL),
51 fHResTPCv0A3(NULL),
52 fHResTPCv0C3(NULL),
53 fHResv0Cv0A3(NULL),
54 fPhiRPv0A(NULL),
55 fPhiRPv0C(NULL),
56 fPhiRPv0Av3(NULL),
57 fPhiRPv0Cv3(NULL),
587d006a 58 fQA(NULL),
59 fQA2(NULL),
60 fQAv3(NULL),
61 fQA2v3(NULL),
afa8df58 62 fPID(new AliFlowBayesianPID()),
587d006a 63 fTree(NULL),
64 fCentrality(-1),
afa8df58 65 evPlAngV0ACor2(0),
66 evPlAngV0CCor2(0),
67 evPlAng2(0),
68 evPlAngV0ACor3(0),
69 evPlAngV0CCor3(0),
70 evPlAng3(0),
587d006a 71 fContAllChargesV0A(NULL),
72 fContAllChargesV0C(NULL),
73 fContAllChargesV0Av3(NULL),
74 fContAllChargesV0Cv3(NULL),
75 fContAllChargesMC(NULL),
243fbce7 76 fContAllChargesMCv3(NULL)
afa8df58 77{
243fbce7 78// DefineOutput(1, TList::Class());
79// DefineOutput(2, TList::Class());
80// DefineOutput(3, TList::Class());
81// DefineOutput(4, TList::Class());
afa8df58 82
83 // Default constructor (should not be used)
84 fList->SetName("resultsV2");
85 fList2->SetName("resultsV3");
587d006a 86 fList3->SetName("resultsMC");
87 fList4->SetName("QA");
afa8df58 88
243fbce7 89 fList->SetOwner(kTRUE);
90 fList2->SetOwner(kTRUE);
91 fList3->SetOwner(kTRUE);
92 fList4->SetOwner(kTRUE);
93
afa8df58 94 fPID->SetNewTrackParam(); // Better tuning for TOF PID tracking effect in LHC10h
afa8df58 95}
96
97//______________________________________________________________________________
98AliAnalysisTaskVnV0::AliAnalysisTaskVnV0(const char *name):
99 AliAnalysisTaskSE(name),
587d006a 100 fVtxCut(10.0), // cut on |vertex| < fVtxCut
afa8df58 101 fEtaCut(0.8), // cut on |eta| < fEtaCut
102 fMinPt(0.15), // cut on pt > fMinPt
243fbce7 103 fV2(kTRUE),
104 fV3(kTRUE),
105 fIsMC(kFALSE),
106 fQAsw(kFALSE),
afa8df58 107 fRun(-1),
108 fList(new TList()),
109 fList2(new TList()),
587d006a 110 fList3(new TList()),
111 fList4(new TList()),
112 fMultV0(NULL),
afa8df58 113 fV0Cpol(100),
114 fV0Apol(100),
587d006a 115 fHResTPCv0A2(NULL),
116 fHResTPCv0C2(NULL),
117 fHResv0Cv0A2(NULL),
118 fHResTPCv0A3(NULL),
119 fHResTPCv0C3(NULL),
120 fHResv0Cv0A3(NULL),
121 fPhiRPv0A(NULL),
122 fPhiRPv0C(NULL),
123 fPhiRPv0Av3(NULL),
124 fPhiRPv0Cv3(NULL),
587d006a 125 fQA(NULL),
126 fQA2(NULL),
127 fQAv3(NULL),
128 fQA2v3(NULL),
afa8df58 129 fPID(new AliFlowBayesianPID()),
587d006a 130 fTree(NULL),
131 fCentrality(-1),
afa8df58 132 evPlAngV0ACor2(0),
133 evPlAngV0CCor2(0),
134 evPlAng2(0),
135 evPlAngV0ACor3(0),
136 evPlAngV0CCor3(0),
137 evPlAng3(0),
587d006a 138 fContAllChargesV0A(NULL),
139 fContAllChargesV0C(NULL),
140 fContAllChargesV0Av3(NULL),
141 fContAllChargesV0Cv3(NULL),
142 fContAllChargesMC(NULL),
243fbce7 143 fContAllChargesMCv3(NULL)
afa8df58 144{
145
146 DefineOutput(1, TList::Class());
147 DefineOutput(2, TList::Class());
587d006a 148 DefineOutput(3, TList::Class());
149 DefineOutput(4, TList::Class());
afa8df58 150
151 // Output slot #1 writes into a TTree
152 fList->SetName("resultsV2");
153 fList2->SetName("resultsV3");
587d006a 154 fList3->SetName("resultsMC");
155 fList4->SetName("QA");
afa8df58 156
243fbce7 157 fList->SetOwner(kTRUE);
158 fList2->SetOwner(kTRUE);
159 fList3->SetOwner(kTRUE);
160 fList4->SetOwner(kTRUE);
161
afa8df58 162 fPID->SetNewTrackParam(); // Better tuning for TOF PID tracking effect in LHC10h
163}
afa8df58 164//_____________________________________________________________________________
165AliAnalysisTaskVnV0::~AliAnalysisTaskVnV0()
166{
167
168}
169
170//______________________________________________________________________________
171void AliAnalysisTaskVnV0::UserCreateOutputObjects()
172{
173
587d006a 174 if(fIsMC) fPID->SetMC(kTRUE);
175
176
afa8df58 177 // Tree for EP debug (comment the adding to v2 list id not needed)
178 fTree = new TTree("tree","tree");
179 fTree->Branch("evPlAngV0ACor2",&evPlAngV0ACor2,"evPlAngV0ACor2/F");
180 fTree->Branch("evPlAngV0CCor2",&evPlAngV0CCor2,"evPlAngV0CCor2/F");
181 fTree->Branch("evPlAng2",&evPlAng2,"evPlAng2/F");
182 fTree->Branch("fCentrality",&fCentrality,"fCentrality/F");
183 fTree->Branch("evPlAngV0ACor3",&evPlAngV0ACor3,"evPlAngV0ACor3/F");
184 fTree->Branch("evPlAngV0CCor3",&evPlAngV0CCor3,"evPlAngV0CCor3/F");
185 fTree->Branch("evPlAng3",&evPlAng3,"evPlAng3/F");
186
187
188 // Container analyses (different steps mean different species)
189 const Int_t nPtBinsTOF = 45;
190 Double_t binsPtTOF[nPtBinsTOF+1] = {0., 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45, 0.5, 0.6, 0.7, 0.8, 0.9, 1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2.0, 2.25, 2.5, 2.75,3.0,3.25,3.5,3.75,4.0,4.5,5,5.5,6,6.5,7,8,9,10,12,15,20};
afa8df58 191 const Int_t nCentrTOF = 9;
243fbce7 192 const Int_t nPsiTOF = 10;
587d006a 193 const Int_t nChargeBinsTOFres = 2;
194 const Int_t nCentrTOFres = 9;
195 const Int_t nProbTOFres = 4;
196 const Int_t nPsiTOFres = 10;
243fbce7 197 const Int_t nMaskPID = 3;
587d006a 198
243fbce7 199 Int_t binsTOF[5] = {nCentrTOFres,nChargeBinsTOFres,nProbTOFres,nPsiTOFres,nMaskPID};
587d006a 200 Int_t binsTOFmc[5] = {nCentrTOFres,nChargeBinsTOFres,1,nPsiTOFres,2};
afa8df58 201
202 // v2 container
587d006a 203 fContAllChargesV0A = new AliFlowVZEROResults("v2A",5,binsTOF);
204 fContAllChargesV0A->SetVarRange(0,-0.5,8.5); // centrality
205 fContAllChargesV0A->SetVarRange(1,-1.5,1.5); // charge
206 fContAllChargesV0A->SetVarRange(2,0.6,1.0001);// prob
207 fContAllChargesV0A->SetVarRange(3,-TMath::Pi()/2,TMath::Pi()/2); // Psi
208 fContAllChargesV0A->SetVarRange(4,-0.5,2.5); // pid mask
209 fContAllChargesV0A->SetVarName(0,"centrality");
210 fContAllChargesV0A->SetVarName(1,"charge");
211 fContAllChargesV0A->SetVarName(2,"prob");
212 fContAllChargesV0A->SetVarName(3,"#Psi");
213 fContAllChargesV0A->SetVarName(4,"PIDmask");
243fbce7 214 if(fV2) fContAllChargesV0A->AddSpecies("all",nPtBinsTOF,binsPtTOF);
215 if(fV2) fContAllChargesV0A->AddSpecies("pi",nPtBinsTOF,binsPtTOF);
216 if(fV2) fContAllChargesV0A->AddSpecies("k",nPtBinsTOF,binsPtTOF);
217 if(fV2) fContAllChargesV0A->AddSpecies("pr",nPtBinsTOF,binsPtTOF);
218 if(fV2) fContAllChargesV0A->AddSpecies("e",nPtBinsTOF,binsPtTOF);
219 if(fV2) fContAllChargesV0A->AddSpecies("d",nPtBinsTOF,binsPtTOF);
220 if(fV2) fContAllChargesV0A->AddSpecies("t",nPtBinsTOF,binsPtTOF);
221 if(fV2) fContAllChargesV0A->AddSpecies("he3",nPtBinsTOF,binsPtTOF);
222 if(fV2) fContAllChargesV0A->AddSpecies("mu",nPtBinsTOF,binsPtTOF);
587d006a 223
224 fContAllChargesV0C = new AliFlowVZEROResults("v2C",5,binsTOF);
225 fContAllChargesV0C->SetVarRange(0,-0.5,8.5); // centrality
226 fContAllChargesV0C->SetVarRange(1,-1.5,1.5); // charge
227 fContAllChargesV0C->SetVarRange(2,0.6,1.0001);// prob
228 fContAllChargesV0C->SetVarRange(3,-TMath::Pi()/2,TMath::Pi()/2); // Psi
229 fContAllChargesV0C->SetVarRange(4,-0.5,2.5); // pid mask
230 fContAllChargesV0C->SetVarName(0,"centrality");
231 fContAllChargesV0C->SetVarName(1,"charge");
232 fContAllChargesV0C->SetVarName(2,"prob");
233 fContAllChargesV0C->SetVarName(3,"#Psi");
234 fContAllChargesV0C->SetVarName(4,"PIDmask");
243fbce7 235 if(fV2) fContAllChargesV0C->AddSpecies("all",nPtBinsTOF,binsPtTOF);
236 if(fV2) fContAllChargesV0C->AddSpecies("pi",nPtBinsTOF,binsPtTOF);
237 if(fV2) fContAllChargesV0C->AddSpecies("k",nPtBinsTOF,binsPtTOF);
238 if(fV2) fContAllChargesV0C->AddSpecies("pr",nPtBinsTOF,binsPtTOF);
239 if(fV2) fContAllChargesV0C->AddSpecies("e",nPtBinsTOF,binsPtTOF);
240 if(fV2) fContAllChargesV0C->AddSpecies("d",nPtBinsTOF,binsPtTOF);
241 if(fV2) fContAllChargesV0C->AddSpecies("t",nPtBinsTOF,binsPtTOF);
242 if(fV2) fContAllChargesV0C->AddSpecies("he3",nPtBinsTOF,binsPtTOF);
243 if(fV2) fContAllChargesV0C->AddSpecies("mu",nPtBinsTOF,binsPtTOF);
afa8df58 244
245 fList->Add(fContAllChargesV0A);
246 fList->Add(fContAllChargesV0C);
247
243fbce7 248 if(fIsMC && fV2){
587d006a 249 fContAllChargesMC = new AliFlowVZEROResults("v2mc",5,binsTOFmc);
250 fContAllChargesMC->SetVarRange(0,-0.5,8.5); // centrality
251 fContAllChargesMC->SetVarRange(1,-1.5,1.5); // charge
252 fContAllChargesMC->SetVarRange(2,0.6,1.0001);// prob
253 fContAllChargesMC->SetVarRange(3,-TMath::Pi()/2,TMath::Pi()/2); // Psi
254 fContAllChargesMC->SetVarRange(4,-0.5,1.5); // pid mask
255 fContAllChargesMC->SetVarName(0,"centrality");
256 fContAllChargesMC->SetVarName(1,"charge");
257 fContAllChargesMC->SetVarName(2,"prob");
258 fContAllChargesMC->SetVarName(3,"#Psi");
259 fContAllChargesMC->SetVarName(4,"PIDmask");
260 fContAllChargesMC->AddSpecies("all",nPtBinsTOF,binsPtTOF);
261 fContAllChargesMC->AddSpecies("pi",nPtBinsTOF,binsPtTOF);
262 fContAllChargesMC->AddSpecies("k",nPtBinsTOF,binsPtTOF);
263 fContAllChargesMC->AddSpecies("pr",nPtBinsTOF,binsPtTOF);
264 fContAllChargesMC->AddSpecies("e",nPtBinsTOF,binsPtTOF);
265 fContAllChargesMC->AddSpecies("mu",nPtBinsTOF,binsPtTOF);
266 fList3->Add(fContAllChargesMC);
267 }
268
afa8df58 269 // v3 container
587d006a 270 fContAllChargesV0Av3 = new AliFlowVZEROResults("v3A",5,binsTOF);
271 fContAllChargesV0Av3->SetVarRange(0,-0.5,8.5); // centrality
272 fContAllChargesV0Av3->SetVarRange(1,-1.5,1.5); // charge
273 fContAllChargesV0Av3->SetVarRange(2,0.6,1.0001);// prob
274 fContAllChargesV0Av3->SetVarRange(3,-TMath::Pi()/3,TMath::Pi()/3); // Psi
275 fContAllChargesV0Av3->SetVarRange(4,-0.5,2.5); // pid mask
276 fContAllChargesV0Av3->SetVarName(0,"centrality");
277 fContAllChargesV0Av3->SetVarName(1,"charge");
278 fContAllChargesV0Av3->SetVarName(2,"prob");
279 fContAllChargesV0Av3->SetVarName(3,"#Psi");
280 fContAllChargesV0Av3->SetVarName(4,"PIDmask");
243fbce7 281 if(fV3) fContAllChargesV0Av3->AddSpecies("all",nPtBinsTOF,binsPtTOF);
282 if(fV3) fContAllChargesV0Av3->AddSpecies("pi",nPtBinsTOF,binsPtTOF);
283 if(fV3) fContAllChargesV0Av3->AddSpecies("k",nPtBinsTOF,binsPtTOF);
284 if(fV3) fContAllChargesV0Av3->AddSpecies("pr",nPtBinsTOF,binsPtTOF);
285 if(fV3) fContAllChargesV0Av3->AddSpecies("e",nPtBinsTOF,binsPtTOF);
286 if(fV3) fContAllChargesV0Av3->AddSpecies("d",nPtBinsTOF,binsPtTOF);
287 if(fV3) fContAllChargesV0Av3->AddSpecies("t",nPtBinsTOF,binsPtTOF);
288 if(fV3) fContAllChargesV0Av3->AddSpecies("he3",nPtBinsTOF,binsPtTOF);
289 if(fV3) fContAllChargesV0Av3->AddSpecies("mu",nPtBinsTOF,binsPtTOF);
587d006a 290
291 fContAllChargesV0Cv3 = new AliFlowVZEROResults("v3C",5,binsTOF);
292 fContAllChargesV0Cv3->SetVarRange(0,-0.5,8.5); // centrality
293 fContAllChargesV0Cv3->SetVarRange(1,-1.5,1.5); // charge
294 fContAllChargesV0Cv3->SetVarRange(2,0.6,1.0001);// prob
295 fContAllChargesV0Cv3->SetVarRange(3,-TMath::Pi()/3,TMath::Pi()/3); // Psi
296 fContAllChargesV0Cv3->SetVarRange(4,-0.5,2.5); // pid mask
297 fContAllChargesV0Cv3->SetVarName(0,"centrality");
298 fContAllChargesV0Cv3->SetVarName(1,"charge");
299 fContAllChargesV0Cv3->SetVarName(2,"prob");
300 fContAllChargesV0Cv3->SetVarName(3,"#Psi");
301 fContAllChargesV0Cv3->SetVarName(4,"PIDmask");
243fbce7 302 if(fV3) fContAllChargesV0Cv3->AddSpecies("all",nPtBinsTOF,binsPtTOF);
303 if(fV3) fContAllChargesV0Cv3->AddSpecies("pi",nPtBinsTOF,binsPtTOF);
304 if(fV3) fContAllChargesV0Cv3->AddSpecies("k",nPtBinsTOF,binsPtTOF);
305 if(fV3) fContAllChargesV0Cv3->AddSpecies("pr",nPtBinsTOF,binsPtTOF);
306 if(fV3) fContAllChargesV0Cv3->AddSpecies("e",nPtBinsTOF,binsPtTOF);
307 if(fV3) fContAllChargesV0Cv3->AddSpecies("d",nPtBinsTOF,binsPtTOF);
308 if(fV3) fContAllChargesV0Cv3->AddSpecies("t",nPtBinsTOF,binsPtTOF);
309 if(fV3) fContAllChargesV0Cv3->AddSpecies("he3",nPtBinsTOF,binsPtTOF);
310 if(fV3) fContAllChargesV0Cv3->AddSpecies("mu",nPtBinsTOF,binsPtTOF);
afa8df58 311
312 fList2->Add(fContAllChargesV0Av3);
313 fList2->Add(fContAllChargesV0Cv3);
314
315 // TProfile for resolutions 3 subevents (V0A, V0C, TPC)
316 // v2
317 fHResTPCv0A2 = new TProfile("hResTPCv0A2","",9,0,9);
318 fHResTPCv0C2 = new TProfile("hResTPCv0C2","",9,0,9);
319 fHResv0Cv0A2 = new TProfile("hResv0Cv0A2","",9,0,9);
320
321 fList->Add(fHResTPCv0A2);
322 fList->Add(fHResTPCv0C2);
323 fList->Add(fHResv0Cv0A2);
324
325 // v3
326 fHResTPCv0A3 = new TProfile("hResTPCv0A3","",9,0,9);
327 fHResTPCv0C3 = new TProfile("hResTPCv0C3","",9,0,9);
328 fHResv0Cv0A3 = new TProfile("hResv0Cv0A3","",9,0,9);
329
330 fList2->Add(fHResTPCv0A3);
331 fList2->Add(fHResTPCv0C3);
332 fList2->Add(fHResv0Cv0A3);
333
334 // V0A and V0C event plane distributions
335 //v2
243fbce7 336 fPhiRPv0A = new TH2F("fPhiRPv0Av2","#phi distribution of EP VZERO-A;centrality;#phi (rad)",9,0,9,nPsiTOF,-TMath::Pi()/2,TMath::Pi()/2);
337 fPhiRPv0C = new TH2F("fPhiRPv0Cv2","#phi distribution of EP VZERO-C;centrality;#phi (rad)",9,0,9,nPsiTOF,-TMath::Pi()/2,TMath::Pi()/2);
afa8df58 338
339 //v3
243fbce7 340 fPhiRPv0Av3 = new TH2F("fPhiRPv0Av3","#phi distribution of EP VZERO-A;centrality;#phi (rad)",9,0,9,nPsiTOF,-TMath::Pi()/3,TMath::Pi()/3);
341 fPhiRPv0Cv3 = new TH2F("fPhiRPv0Cv3","#phi distribution of EP VZERO-C;centrality;#phi (rad)",9,0,9,nPsiTOF,-TMath::Pi()/3,TMath::Pi()/3);
afa8df58 342
343 // QA container
344 // v2
345 const Int_t nDETsignal = 50;
346 Double_t binDETsignal[nDETsignal+1];
347 for(Int_t i=0;i<nDETsignal+1;i++){
348 binDETsignal[i] = -5 + i*10. / nDETsignal;
349 }
587d006a 350// const Int_t nEta = 5;
351// Double_t binEta[nEta+1];
352// for(Int_t i=0;i<nEta+1;i++){
353// binEta[i] = -1 + i*2. / nEta;
354// }
afa8df58 355
356 const Int_t nDeltaPhi = 5;
afa8df58 357 const Int_t nDeltaPhiV3 = 7;
afa8df58 358
243fbce7 359 Int_t binsQA[5] = {nCentrTOF,7,5,nDeltaPhi,2};
360 Int_t binsQAv3[5] = {nCentrTOF,7,5,nDeltaPhiV3,2};
361
362
363 fQA = new AliFlowVZEROQA("v2AQA",5,binsQA);
364 fQA->SetVarRange(0,-0.5,8.5); // centrality
365 fQA->SetVarRange(1,0,7); // pt
366 fQA->SetVarRange(2,0.,1.0001);// prob
367 fQA->SetVarRange(3,-TMath::Pi(),TMath::Pi()); // Psi
368 fQA->SetVarRange(4,-0.5,1.5); // pid mask
369 fQA->SetVarName(0,"centrality");
370 fQA->SetVarName(1,"p_{t}");
371 fQA->SetVarName(2,"prob");
372 fQA->SetVarName(3,"#Psi");
373 fQA->SetVarName(4,"PIDmask");
374 if(fQAsw && fV2) fQA->AddSpecies("pi",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
375 if(fQAsw && fV2) fQA->AddSpecies("k",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
376 if(fQAsw && fV2) fQA->AddSpecies("pr",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
377// if(fQAsw && fV2) fQA->AddSpecies("e",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
378// fQA->AddSpecies("d",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
379// fQA->AddSpecies("t",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
380// fQA->AddSpecies("he3",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
381
382 fQA2 = new AliFlowVZEROQA("v2CQA",5,binsQA);
383 fQA2->SetVarRange(0,-0.5,8.5); // centrality
384 fQA2->SetVarRange(1,0,7); // pt
385 fQA2->SetVarRange(2,0.,1.0001);// prob
386 fQA2->SetVarRange(3,-TMath::Pi(),TMath::Pi()); // Psi
387 fQA2->SetVarRange(4,-0.5,1.5); // pid mask
388 fQA2->SetVarName(0,"centrality");
389 fQA2->SetVarName(1,"p_{t}");
390 fQA2->SetVarName(2,"prob");
391 fQA2->SetVarName(3,"#Psi");
392 fQA2->SetVarName(4,"PIDmask");
393 if(fQAsw && fV2) fQA2->AddSpecies("pi",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
394 if(fQAsw && fV2) fQA2->AddSpecies("k",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
395 if(fQAsw && fV2) fQA2->AddSpecies("pr",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
396// if(fQAsw && fV2) fQA2->AddSpecies("e",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
397// fQA2->AddSpecies("d",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
398// fQA2->AddSpecies("t",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
399// fQA2->AddSpecies("he3",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
400
401 fQAv3 = new AliFlowVZEROQA("v3AQA",5,binsQAv3);
402 fQAv3->SetVarRange(0,-0.5,8.5); // centrality
403 fQAv3->SetVarRange(1,0,7); // pt
404 fQAv3->SetVarRange(2,0.,1.0001);// prob
405 fQAv3->SetVarRange(3,-TMath::Pi(),TMath::Pi()); // Psi
406 fQAv3->SetVarRange(4,-0.5,1.5); // pid mask
407 fQAv3->SetVarName(0,"centrality");
408 fQAv3->SetVarName(1,"p_{t}");
409 fQAv3->SetVarName(2,"prob");
410 fQAv3->SetVarName(3,"#Psi");
411 fQAv3->SetVarName(4,"PIDmask");
412 if(fQAsw && fV3) fQAv3->AddSpecies("pi",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
413// if(fQAsw && fV3) fQAv3->AddSpecies("k",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
414// if(fQAsw && fV3) fQAv3->AddSpecies("pr",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
415// if(fQAsw && fV2) fQAv3->AddSpecies("e",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
416// fQAv3->AddSpecies("d",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
417// fQAv3->AddSpecies("t",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
418// fQAv3->AddSpecies("he3",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
419
420 fQA2v3 = new AliFlowVZEROQA("v3CQA",5,binsQAv3);
421 fQA2v3->SetVarRange(0,-0.5,8.5); // centrality
422 fQA2v3->SetVarRange(1,0,7); // pt
423 fQA2v3->SetVarRange(2,0.,1.0001);// prob
424 fQA2v3->SetVarRange(3,-TMath::Pi(),TMath::Pi()); // Psi
425 fQA2v3->SetVarRange(4,-0.5,1.5); // pid mask
426 fQA2v3->SetVarName(0,"centrality");
427 fQA2v3->SetVarName(1,"p_{t}");
428 fQA2v3->SetVarName(2,"prob");
429 fQA2v3->SetVarName(3,"#Psi");
430 fQA2v3->SetVarName(4,"PIDmask");
431 if(fQAsw && fV3) fQA2v3->AddSpecies("pi",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
432// if(fQAsw && fV3) fQA2v3->AddSpecies("k",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
433// if(fQAsw && fV3) fQA2v3->AddSpecies("pr",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
434// if(fQAsw && fV2) fQA2v3->AddSpecies("e",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
435// fQA2v3->AddSpecies("d",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
436// fQA2v3->AddSpecies("t",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
437// fQA2v3->AddSpecies("he3",nDETsignal,binDETsignal,nDETsignal,binDETsignal);
afa8df58 438
439 fList->Add(fPhiRPv0A);
440 fList->Add(fPhiRPv0C);
587d006a 441
442 if(fQAsw && fV2){
443 fList4->Add(fQA);
444 fList4->Add(fQA2);
445 }
afa8df58 446
447 fList2->Add(fPhiRPv0Av3);
448 fList2->Add(fPhiRPv0Cv3);
afa8df58 449
587d006a 450 if(fQAsw && fV3){
451 fList4->Add(fQAv3);
452 fList4->Add(fQA2v3);
453 }
454
a72e0394 455 // fList->Add(fTree); // comment if not needed
afa8df58 456
457 printf("Output creation ok!!\n\n\n\n");
458
459 // Post output data.
460 if(fV2) PostData(1, fList);
461 if(fV3) PostData(2, fList2);
587d006a 462 if(fIsMC) PostData(3, fList3);
463 if(fQAsw) PostData(4, fList4);
243fbce7 464
afa8df58 465}
466
467//______________________________________________________________________________
468void AliAnalysisTaskVnV0::UserExec(Option_t *)
469{
470 // Main loop
471 // Called for each event
472
243fbce7 473 fOutputAOD = dynamic_cast<AliAODEvent*>(InputEvent());
474 if(!fOutputAOD){
afa8df58 475 Printf("%s:%d AODEvent not found in Input Manager",(char*)__FILE__,__LINE__);
476 this->Dump();
477 return;
478 }
479
243fbce7 480 Int_t run = fOutputAOD->GetRunNumber();
afa8df58 481
482 if(run != fRun){
483 // Load the calibrations run dependent
484 OpenInfoCalbration(run);
485 fRun=run;
486 }
487
243fbce7 488 Float_t zvtx = GetVertex(fOutputAOD);
afa8df58 489
587d006a 490
491
492 //Get the MC object
d6a1c304 493 if(fIsMC){
243fbce7 494 AliAODMCHeader *mcHeader = dynamic_cast<AliAODMCHeader*>(fOutputAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
d6a1c304 495 if (!mcHeader) {
496 AliError("Could not find MC Header in AOD");
497 return;
498 }
587d006a 499 }
500
501 /*
502 AliMCEvent* mcEvent = MCEvent();
503 if (!mcEvent) {
504 Printf("ERROR: Could not retrieve MC event");
505 return;
506 }
507
508 Double_t gReactionPlane = -999., gImpactParameter = -999.;
509 //Get the MC header
510 AliGenHijingEventHeader* headerH = dynamic_cast<AliGenHijingEventHeader*>(mcEvent->GenEventHeader());
511 if (headerH) {
512 //Printf("=====================================================");
513 //Printf("Reaction plane angle: %lf",headerH->ReactionPlaneAngle());
514 //Printf("=====================================================");
515 gReactionPlane = headerH->ReactionPlaneAngle();
516 gImpactParameter = headerH->ImpactParameter();
517 }
518
519*/
520
afa8df58 521 if (TMath::Abs(zvtx) < fVtxCut) {
522 //Centrality
523 Float_t v0Centr = -10.;
524 Float_t trkCentr = -10.;
243fbce7 525 AliCentrality *centrality = fOutputAOD->GetCentrality();
afa8df58 526 if (centrality){
527 v0Centr = centrality->GetCentralityPercentile("V0M");
528 trkCentr = centrality->GetCentralityPercentile("TRK");
529 }
530
531 if(TMath::Abs(v0Centr - trkCentr) < 5.0){ // consistency cut on centrality selection
243fbce7 532 fPID->SetDetResponse(fOutputAOD, v0Centr); // Set the PID object for each event!!!!
533 Analyze(fOutputAOD,v0Centr); // Do analysis!!!!
afa8df58 534
535 fCentrality = v0Centr;
536 if(fV2) fTree->Fill();
537 }
538 }
539
540}
541
542//________________________________________________________________________
543void AliAnalysisTaskVnV0::Analyze(AliAODEvent* aodEvent, Float_t v0Centr)
544{
545 Float_t mass[8] = {5.10998909999999971e-04, 1.05658000000000002e-01, 1.39570000000000000e-01, 4.93676999999999977e-01, 9.38271999999999995e-01,1.87783699999999998,2.81740199999999996,1.40805449999999999};
546
547 // Event plane resolution for v2
548 Float_t evPlRes[18] = {0.350582,0.505393,0.607845,0.632913,0.592230,0.502489,0.381717,0.249539,0.133180, // V0A vs. centrality
549 0.446480,0.612705,0.712222,0.736200,0.697907,0.610114,0.481009,0.327402,0.182277};// V0C vs. centrality
550
551 Int_t iC = -1;
552 if (v0Centr >0 && v0Centr < 80){ // analysis only for 0-80% centrality classes
553 // centrality bins
554 if(v0Centr < 5) iC = 0;
555 else if(v0Centr < 10) iC = 1;
556 else if(v0Centr < 20) iC = 2;
557 else if(v0Centr < 30) iC = 3;
558 else if(v0Centr < 40) iC = 4;
559 else if(v0Centr < 50) iC = 5;
560 else if(v0Centr < 60) iC = 6;
561 else if(v0Centr < 70) iC = 7;
562 else iC = 8;
563
564 //reset Q vector info
565 Double_t Qxa2 = 0, Qya2 = 0;
566 Double_t Qxc2 = 0, Qyc2 = 0;
567 Double_t Qxa3 = 0, Qya3 = 0;
568 Double_t Qxc3 = 0, Qyc3 = 0;
587d006a 569
570 Int_t nAODTracks = aodEvent->GetNumberOfTracks();
571
572 AliAODMCHeader *mcHeader = NULL;
573 TClonesArray *mcArray = NULL;
574 Float_t evplaneMC = 0;
575 if(fIsMC){
243fbce7 576 mcHeader = dynamic_cast<AliAODMCHeader*>(fOutputAOD->GetList()->FindObject(AliAODMCHeader::StdBranchName()));
587d006a 577
578 if (mcHeader) {
4edd7307 579 evplaneMC = mcHeader->GetReactionPlaneAngle();
580 if(evplaneMC > TMath::Pi()/2 && evplaneMC <= TMath::Pi()*3/2) evplaneMC-=TMath::Pi();
581 else if(evplaneMC > TMath::Pi()*3/2) evplaneMC-=2*TMath::Pi();
243fbce7 582 mcArray = (TClonesArray*)fOutputAOD->GetList()->FindObject(AliAODMCParticle::StdBranchName());
587d006a 583 }
584 }
585
afa8df58 586 //V0 info
587 AliAODVZERO* aodV0 = aodEvent->GetVZEROData();
588
589 for (Int_t iv0 = 0; iv0 < 64; iv0++) {
590 Double_t phiV0 = TMath::PiOver4()*(0.5 + iv0 % 8);
591 Float_t multv0 = aodV0->GetMultiplicity(iv0);
592 if (iv0 < 32){ // V0C
593 Qxc2 += TMath::Cos(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
594 Qyc2 += TMath::Sin(2*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
595 Qxc3 += TMath::Cos(3*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
596 Qyc3 += TMath::Sin(3*phiV0) * multv0*fV0Cpol/fMultV0->GetBinContent(iv0+1);
597 } else { // V0A
598 Qxa2 += TMath::Cos(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
599 Qya2 += TMath::Sin(2*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
600 Qxa3 += TMath::Cos(3*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
601 Qya3 += TMath::Sin(3*phiV0) * multv0*fV0Apol/fMultV0->GetBinContent(iv0+1);
602 }
603 }
604
605 //grab for each centrality the proper histo with the Qx and Qy to do the recentering
606 Double_t Qxamean2 = fMeanQ[iC][1][0];
607 Double_t Qxarms2 = fWidthQ[iC][1][0];
608 Double_t Qyamean2 = fMeanQ[iC][1][1];
609 Double_t Qyarms2 = fWidthQ[iC][1][1];
610 Double_t Qxamean3 = fMeanQv3[iC][1][0];
611 Double_t Qxarms3 = fWidthQv3[iC][1][0];
612 Double_t Qyamean3 = fMeanQv3[iC][1][1];
613 Double_t Qyarms3 = fWidthQv3[iC][1][1];
614
615 Double_t Qxcmean2 = fMeanQ[iC][0][0];
616 Double_t Qxcrms2 = fWidthQ[iC][0][0];
617 Double_t Qycmean2 = fMeanQ[iC][0][1];
618 Double_t Qycrms2 = fWidthQ[iC][0][1];
619 Double_t Qxcmean3 = fMeanQv3[iC][0][0];
620 Double_t Qxcrms3 = fWidthQv3[iC][0][0];
621 Double_t Qycmean3 = fMeanQv3[iC][0][1];
622 Double_t Qycrms3 = fWidthQv3[iC][0][1];
623
624 Double_t QxaCor2 = (Qxa2 - Qxamean2)/Qxarms2;
625 Double_t QyaCor2 = (Qya2 - Qyamean2)/Qyarms2;
626 Double_t QxcCor2 = (Qxc2 - Qxcmean2)/Qxcrms2;
627 Double_t QycCor2 = (Qyc2 - Qycmean2)/Qycrms2;
628 Double_t QxaCor3 = (Qxa3 - Qxamean3)/Qxarms3;
629 Double_t QyaCor3 = (Qya3 - Qyamean3)/Qyarms3;
630 Double_t QxcCor3 = (Qxc3 - Qxcmean3)/Qxcrms3;
631 Double_t QycCor3 = (Qyc3 - Qycmean3)/Qycrms3;
632
633 evPlAngV0ACor2 = TMath::ATan2(QyaCor2, QxaCor2)/2.;
634 evPlAngV0CCor2 = TMath::ATan2(QycCor2, QxcCor2)/2.;
635 evPlAngV0ACor3 = TMath::ATan2(QyaCor3, QxaCor3)/3.;
636 evPlAngV0CCor3 = TMath::ATan2(QycCor3, QxcCor3)/3.;
637
afa8df58 638 //loop track and get pid
639 for(Int_t iT = 0; iT < nAODTracks; iT++) { // loop on the tracks
640 AliAODTrack* aodTrack = aodEvent->GetTrack(iT);
641
642 if (!aodTrack){
643 aodTrack->Delete();
644 continue;
645 }
646
647 Bool_t trkFlag = aodTrack->TestFilterBit(1); // TPC only tracks
648
649 if ((TMath::Abs(aodTrack->Eta()) > fEtaCut) || (aodTrack->Pt() < fMinPt) || (aodTrack->GetTPCNcls() < 70) || !trkFlag){
650 continue;
651 }
652
653 Double_t b[2] = {-99., -99.};
654 Double_t bCov[3] = {-99., -99., -99.};
243fbce7 655 if (!aodTrack->PropagateToDCA(fOutputAOD->GetPrimaryVertex(), fOutputAOD->GetMagneticField(), 100., b, bCov))
afa8df58 656 continue;
657
658 if ((TMath::Abs(b[0]) > 3.0) || (TMath::Abs(b[1]) > 2.4))
659 continue;
660
661 // re-map the container in an array to do the analysis for V0A and V0C within a loop
662 Float_t evPlAngV0[2] = {evPlAngV0ACor2,evPlAngV0CCor2};
587d006a 663 AliFlowVZEROResults *contV0[2] = {fContAllChargesV0A,fContAllChargesV0C};
243fbce7 664 AliFlowVZEROQA *QA[2] = {fQA,fQA2};
afa8df58 665
666 Float_t evPlAngV0v3[2] = {evPlAngV0ACor3,evPlAngV0CCor3};
587d006a 667 AliFlowVZEROResults *contV0v3[2] = {fContAllChargesV0Av3,fContAllChargesV0Cv3};
243fbce7 668 AliFlowVZEROQA *QAv3[2] = {fQAv3,fQA2v3};
afa8df58 669
587d006a 670 // Fill MC results
671 if(fIsMC && mcArray){
243fbce7 672 fPID->ComputeProb(aodTrack,fOutputAOD); // compute Bayesian probabilities
587d006a 673 Float_t tofMismProbMC = fPID->GetTOFMismProb(); // TOF mismatch probability requested to be lower than 50% for TOF analysis
674
675 Float_t xMC[5] = {iC,aodTrack->Charge(),1,evplaneMC,fPID->GetCurrentMask(1)&&tofMismProbMC < 0.5}; // to fill analysis v2 container
676
677
678 Float_t v2mc = TMath::Cos(2*(aodTrack->Phi() - evplaneMC));
679
680 fContAllChargesMC->Fill(0,aodTrack->Pt(),v2mc,xMC);
681
682 Int_t iS = TMath::Abs(((AliAODMCParticle*)mcArray->At(TMath::Abs(aodTrack->GetLabel())))->GetPdgCode());
683 if(iS==11){
684 fContAllChargesMC->Fill(4,aodTrack->Pt(),v2mc,xMC);
685 }
686 else if(iS==13){
687 fContAllChargesMC->Fill(5,aodTrack->Pt(),v2mc,xMC);
688 }
689 else if(iS==211){
690 fContAllChargesMC->Fill(1,aodTrack->Pt(),v2mc,xMC);
691 }
692 else if(iS==321){
693 fContAllChargesMC->Fill(2,aodTrack->Pt(),v2mc,xMC);
694 }
695 else if(iS==2212){
696 fContAllChargesMC->Fill(3,aodTrack->Pt(),v2mc,xMC);
697 }
698 }
699
afa8df58 700 for(Int_t iV0=0;iV0<2;iV0++){ // loop on A and C side
701
702 fPID->SetPsiCorrectionDeDx(evPlAngV0[iV0],evPlRes[iV0*8+iC]); // set the PID dE/dx correction as a function of the v2-EP (resolution is needed)
703
587d006a 704 Float_t v2V0 = TMath::Cos(2*(aodTrack->Phi() - evPlAngV0[iV0]));
705 Float_t v3V0 = TMath::Cos(3*(aodTrack->Phi() - evPlAngV0v3[iV0]));
afa8df58 706
243fbce7 707 fPID->ComputeProb(aodTrack,fOutputAOD); // compute Bayesian probabilities
587d006a 708 Float_t dedx = fPID->GetDeDx();//aodTrack->GetTPCsignal();
afa8df58 709 Float_t *probRead = fPID->GetProb();
710 Float_t prob[8] = {probRead[0],probRead[1],probRead[2],probRead[3],probRead[4],probRead[5],probRead[6],probRead[7]};
711 Float_t tofMismProb = fPID->GetTOFMismProb(); // TOF mismatch probability requested to be lower than 50% for TOF analysis
587d006a 712 Float_t x[5] = {iC,aodTrack->Charge(),1,evPlAngV0[iV0],fPID->GetCurrentMask(1)&&tofMismProb < 0.5}; // to fill analysis v2 container
713 Float_t x3[5] = {iC,aodTrack->Charge(),1,evPlAngV0v3[iV0],fPID->GetCurrentMask(1)&&tofMismProb < 0.5}; // to fill analysis v3 container
afa8df58 714
afa8df58 715 // Fill no PID
587d006a 716 if(fV2) contV0[iV0]->Fill(0,aodTrack->Pt(),v2V0,x);
717 if(fV3) contV0v3[iV0]->Fill(0,aodTrack->Pt(),v3V0,x3);
718
afa8df58 719
afa8df58 720 Double_t dedxExp[8];
721 Float_t tof = -1;
722 Double_t inttimes[8] = {0.,0.,0.,0.,0.,0.,0.,0.};
723 Double_t expTOFsigma[8] = {0.,0.,0.,0.,0.,0.,0.,0.};
78a615f6 724
725 Float_t nsigmaTPC[8];
726 Float_t nsigmaTOF[8];
727
afa8df58 728 if(aodTrack->GetDetPid()){ // check the PID object is available
78a615f6 729 for(Int_t iS=0;iS < 8;iS++){
afa8df58 730 dedxExp[iS] = fPID->GetExpDeDx(aodTrack,iS);
78a615f6 731 nsigmaTPC[iS] = (dedx - dedxExp[iS])/(dedxExp[iS]*0.07);
732 // printf("TPC %i = %f (%f %f)\n",iS, nsigmaTPC[iS],dedx, dedxExp[iS]);
733 }
587d006a 734
afa8df58 735 if(fPID->GetCurrentMask(1)){ // if TOF is present
736 Float_t ptrack = aodTrack->P();
737 tof = aodTrack->GetTOFsignal() - fPID->GetESDpid()->GetTOFResponse().GetStartTime(ptrack);
738 aodTrack->GetIntegratedTimes(inttimes);
739
740 for(Int_t iS=5;iS < 8;iS++) // extra info for light nuclei
741 inttimes[iS] = inttimes[0] / ptrack * mass[iS] * TMath::Sqrt(1+ptrack*ptrack/mass[iS]/mass[iS]);
742
78a615f6 743 for(Int_t iS=0;iS<8;iS++){
744 expTOFsigma[iS] = fPID->GetESDpid()->GetTOFResponse().GetExpectedSigma(ptrack, inttimes[iS], mass[iS]);
745 nsigmaTOF[iS] = (tof - inttimes[iS])/expTOFsigma[iS];
746 // printf("TOF %i = %f\n",iS, nsigmaTOF[iS]);
747 }
afa8df58 748 }
749 }
750
751 Float_t deltaPhiV0 = aodTrack->Phi() - evPlAngV0[iV0];
752 if(deltaPhiV0 > TMath::Pi()) deltaPhiV0 -= 2*TMath::Pi();
753 else if(deltaPhiV0 < -TMath::Pi()) deltaPhiV0 += 2*TMath::Pi();
754 if(deltaPhiV0 > TMath::Pi()) deltaPhiV0 -= 2*TMath::Pi();
755 else if(deltaPhiV0 < -TMath::Pi()) deltaPhiV0 += 2*TMath::Pi();
756
757 Float_t deltaPhiV0v3 = aodTrack->Phi() - evPlAngV0v3[iV0];
758 if(deltaPhiV0v3 > TMath::Pi()) deltaPhiV0v3 -= 2*TMath::Pi();
759 else if(deltaPhiV0v3 < -TMath::Pi()) deltaPhiV0v3 += 2*TMath::Pi();
760 if(deltaPhiV0v3 > TMath::Pi()) deltaPhiV0v3 -= 2*TMath::Pi();
761 else if(deltaPhiV0v3 < -TMath::Pi()) deltaPhiV0v3 += 2*TMath::Pi();
762
763 // variable to fill QA container
243fbce7 764 Float_t xQA[5] = {iC,aodTrack->Pt(), 0.0,deltaPhiV0,x[4]}; // v2
765 Float_t xQA3[5] = {iC,aodTrack->Pt(), 0.0,deltaPhiV0v3,x[4]}; // v3
afa8df58 766
78a615f6 767 // QA fill
243fbce7 768 if(!(fPID->GetCurrentMask(0)) || !aodTrack->GetDetPid() || dedx < 10. || aodTrack->Pt() < 0 || aodTrack->Pt() > 7){}
769 else{
770 if(TMath::Abs(nsigmaTPC[2])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[2])<5))){ //pi
771 xQA[2] = prob[2];
772 xQA3[2] = xQA[2];
773 if(fV2) QA[iV0]->Fill(0,nsigmaTPC[2],nsigmaTOF[2],xQA);
774 if(fV3) QAv3[iV0]->Fill(0,nsigmaTPC[2],nsigmaTOF[2],xQA3);
78a615f6 775 }
243fbce7 776 if(TMath::Abs(nsigmaTPC[3])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[3])<5))){ //K
777 xQA[2] = prob[3];
778 xQA3[2] = xQA[2];
779 if(fV2) QA[iV0]->Fill(1,nsigmaTPC[3],nsigmaTOF[3],xQA);
780// if(fV3) QAv3[iV0]->Fill(1,nsigmaTPC[3],nsigmaTOF[3],xQA3);
78a615f6 781 }
243fbce7 782 if(TMath::Abs(nsigmaTPC[4])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[4])<5))){//p
783 xQA[2] = prob[4];
784 xQA3[2] = xQA[2];
785 if(fV2) QA[iV0]->Fill(2,nsigmaTPC[4],nsigmaTOF[4],xQA);
786// if(fV3) QAv3[iV0]->Fill(2,nsigmaTPC[4],nsigmaTOF[4],xQA3);
78a615f6 787 }
243fbce7 788 if(TMath::Abs(nsigmaTPC[0])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[0])<5))){//e
789 xQA[2] = prob[0];
790 xQA3[2] = xQA[2];
791// if(fV2) QA[iV0]->Fill(3,nsigmaTPC[0],nsigmaTOF[0],xQA);
792// if(fV3) QAv3[iV0]->Fill(3,nsigmaTPC[0],nsigmaTOF[0],xQA3);
78a615f6 793 }
243fbce7 794 if(TMath::Abs(nsigmaTPC[5])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[5])<5))){//d
795 xQA[2] = prob[5];
796 xQA3[2] = xQA[2];
797 // if(fV2) QA[iV0]->Fill(4,nsigmaTPC[5],nsigmaTOF[5],xQA);
798 // if(fV3) QAv3[iV0]->Fill(4,nsigmaTPC[5],nsigmaTOF[5],xQA3);
78a615f6 799 }
243fbce7 800 if(TMath::Abs(nsigmaTPC[6])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[6])<5))){//t
801 xQA[2] = prob[6];
802 xQA3[2] = xQA[2];
803 // if(fV2) QA[iV0]->Fill(5,nsigmaTPC[6],nsigmaTOF[6],xQA);
804 // if(fV3) QAv3[iV0]->Fill(5,nsigmaTPC[6],nsigmaTOF[6],xQA3);
78a615f6 805 }
243fbce7 806 if(TMath::Abs(nsigmaTPC[7])<5 && (!(fPID->GetCurrentMask(1)) || (TMath::Abs(nsigmaTOF[7])<5))){//He3
807 xQA[2] = prob[7];
808 xQA3[2] = xQA[2];
809 // if(fV2) QA[iV0]->Fill(6,nsigmaTPC[7],nsigmaTOF[7],xQA);
810 // if(fV3) QAv3[iV0]->Fill(6,nsigmaTPC[7],nsigmaTOF[7],xQA3);
78a615f6 811 }
78a615f6 812 }
813
afa8df58 814 //pid selection
815 if(!(fPID->GetCurrentMask(0)) || !aodTrack->GetDetPid()){} // TPC PID and PID object strictly required (very important!!!!)
816 else if(prob[2] > 0.6){ // pi
587d006a 817 x[2] = prob[2];
587d006a 818 x3[2] = x[2];
78a615f6 819 if(TMath::Abs(nsigmaTPC[2]) < 5){ // TPC 5 sigma extra cut to accept the track
587d006a 820 if(fV2) contV0[iV0]->Fill(1,aodTrack->Pt(),v2V0,x);
a72e0394 821 if(fV3) contV0v3[iV0]->Fill(1,aodTrack->Pt(),v3V0,x3);
afa8df58 822 }
823 }
824 else if(prob[3] > 0.6){ // K
587d006a 825 x[2] = prob[3];
587d006a 826 x3[2] = x[2];
78a615f6 827 if(TMath::Abs(nsigmaTPC[3]) < 5){
587d006a 828 if(fV2) contV0[iV0]->Fill(2,aodTrack->Pt(),v2V0,x);
a72e0394 829 if(fV3) contV0v3[iV0]->Fill(2,aodTrack->Pt(),v3V0,x3);
afa8df58 830 }
831 }
832 else if(prob[4] > 0.6){ // p
587d006a 833 x[2] = prob[4];
587d006a 834 x3[2] = x[2];
78a615f6 835 if(TMath::Abs(nsigmaTPC[4]) < 5){
587d006a 836 if(fV2) contV0[iV0]->Fill(3,aodTrack->Pt(),v2V0,x);
a72e0394 837 if(fV3) contV0v3[iV0]->Fill(3,aodTrack->Pt(),v3V0,x3);
afa8df58 838 }
839 }
840 else if(prob[0] > 0.6){ // e
587d006a 841 x[2] = prob[0];
587d006a 842 x3[2] = x[2];
78a615f6 843 if(TMath::Abs(nsigmaTPC[0]) < 5){
587d006a 844 if(fV2) contV0[iV0]->Fill(4,aodTrack->Pt(),v2V0,x);
a72e0394 845 if(fV3) contV0v3[iV0]->Fill(4,aodTrack->Pt(),v3V0,x3);
afa8df58 846 }
847 }
587d006a 848 else if(prob[1] > 0.6){ // mu
587d006a 849 x[2] = prob[1];
587d006a 850 x3[2] = x[2];
78a615f6 851 if(TMath::Abs(nsigmaTPC[1]) < 5){
587d006a 852 if(fV2) contV0[iV0]->Fill(8,aodTrack->Pt(),v2V0,x);
a72e0394 853 if(fV3) contV0v3[iV0]->Fill(8,aodTrack->Pt(),v3V0,x3);
587d006a 854 }
855 }
afa8df58 856 else if(prob[5] > 0.6){ // d
587d006a 857 x[2] = prob[5];
587d006a 858 x3[2] = x[2];
78a615f6 859 if(TMath::Abs(nsigmaTPC[5]) < 5){
587d006a 860 if(fV2) contV0[iV0]->Fill(5,aodTrack->Pt(),v2V0,x);
a72e0394 861 if(fV3) contV0v3[iV0]->Fill(5,aodTrack->Pt(),v3V0,x3);
afa8df58 862 }
863 }
864 else if(prob[6] > 0.6){ // t
587d006a 865 x[2] = prob[6];
587d006a 866 x3[2] = x[2];
78a615f6 867 if(TMath::Abs(nsigmaTPC[6]) < 5){
587d006a 868 if(fV2) contV0[iV0]->Fill(6,aodTrack->Pt(),v2V0,x);
a72e0394 869 if(fV3) contV0v3[iV0]->Fill(6,aodTrack->Pt(),v3V0,x3);
afa8df58 870 }
871 }
872 else if(prob[7] > 0.6){ // He3
587d006a 873 x[2] = prob[7];
587d006a 874 x3[2] = x[2];
78a615f6 875 if(TMath::Abs(nsigmaTPC[7]) < 5){
587d006a 876 if(fV2) contV0[iV0]->Fill(7,aodTrack->Pt()*2,v2V0,x);
a72e0394 877 if(fV3) contV0v3[iV0]->Fill(7,aodTrack->Pt()*2,v3V0,x3);
afa8df58 878 }
afa8df58 879 }
880
587d006a 881 if(x[4] > 0.5){ // if TOF was present redo TPC stand alone PID to check the PID in the same acceptance (PID mask = 2)
afa8df58 882 fPID->ResetDetOR(1); // exclude TOF from PID
883 tofMismProb = 0;
884
243fbce7 885 fPID->ComputeProb(aodTrack,fOutputAOD);
587d006a 886 dedx = fPID->GetDeDx();//aodTrack->GetTPCsignal();
afa8df58 887 probRead = fPID->GetProb();
888
889 fPID->SetDetOR(1); // include TOF for PID
890 }
891 Float_t probTPC[8] = {probRead[0],probRead[1],probRead[2],probRead[3],probRead[4],probRead[5],probRead[6],probRead[7]}; // TPC stand alone prbabilities
892
893 //pid selection TPC S.A. with TOF matching
587d006a 894 x[4]*=2; // set the mask to 2 id TOF is present
895 if(x[4]<1 || !(fPID->GetCurrentMask(0)) || !aodTrack->GetDetPid()){} // TPC PID S.A. PID in TOF acceptance
afa8df58 896 else if(probTPC[2] > 0.6){ // pi
587d006a 897 x[2] = probTPC[2];
898 x3[2] = x[2];
78a615f6 899 if(TMath::Abs(nsigmaTPC[2]) < 5){
587d006a 900 if(fV2) contV0[iV0]->Fill(1,aodTrack->Pt(),v2V0,x);
a72e0394 901 if(fV3) contV0v3[iV0]->Fill(1,aodTrack->Pt(),v3V0,x3);
afa8df58 902 }
903 }
904 else if(probTPC[3] > 0.6){ // K
587d006a 905 x[2] = probTPC[3];
906 x3[2] = x[2];
78a615f6 907 if(TMath::Abs(nsigmaTPC[3]) < 5){
587d006a 908 if(fV2) contV0[iV0]->Fill(2,aodTrack->Pt(),v2V0,x);
a72e0394 909 if(fV3) contV0v3[iV0]->Fill(2,aodTrack->Pt(),v3V0,x3);
afa8df58 910 }
911 }
912 else if(probTPC[4] > 0.6){ // p
587d006a 913 x[2] = probTPC[4];
914 x3[2] = x[2];
78a615f6 915 if(TMath::Abs(nsigmaTPC[4]) < 5){
587d006a 916 if(fV2) contV0[iV0]->Fill(3,aodTrack->Pt(),v2V0,x);
a72e0394 917 if(fV3) contV0v3[iV0]->Fill(3,aodTrack->Pt(),v3V0,x3);
afa8df58 918 }
919 }
920 else if(probTPC[0] > 0.6){ // e
587d006a 921 x[2] = probTPC[0];
922 x3[2] = x[2];
78a615f6 923 if(TMath::Abs(nsigmaTPC[0]) < 5){
587d006a 924 if(fV2) contV0[iV0]->Fill(4,aodTrack->Pt(),v2V0,x);
a72e0394 925 if(fV3) contV0v3[iV0]->Fill(4,aodTrack->Pt(),v3V0,x3);
587d006a 926 }
927 }
928 else if(probTPC[1] > 0.6){ // mu
929 x[2] = probTPC[1];
930 x3[2] = x[2];
78a615f6 931 if(TMath::Abs(nsigmaTPC[1]) < 5){
587d006a 932 if(fV2) contV0[iV0]->Fill(8,aodTrack->Pt(),v2V0,x);
a72e0394 933 if(fV3) contV0v3[iV0]->Fill(8,aodTrack->Pt(),v3V0,x3);
afa8df58 934 }
935 }
936 else if(probTPC[5] > 0.6){ // d
587d006a 937 x[2] = probTPC[5];
938 x3[2] = x[2];
78a615f6 939 if(TMath::Abs(nsigmaTPC[5]) < 5){
587d006a 940 if(fV2) contV0[iV0]->Fill(5,aodTrack->Pt(),v2V0,x);
a72e0394 941 if(fV3) contV0v3[iV0]->Fill(5,aodTrack->Pt(),v3V0,x3);
afa8df58 942 }
943 }
944 else if(probTPC[6] > 0.6){ // t
587d006a 945 x[2] = probTPC[6];
946 x3[2] = x[2];
78a615f6 947 if(TMath::Abs(nsigmaTPC[6]) < 5){
587d006a 948 if(fV2) contV0[iV0]->Fill(6,aodTrack->Pt(),v2V0,x);
a72e0394 949 if(fV3) contV0v3[iV0]->Fill(6,aodTrack->Pt(),v3V0,x3);
afa8df58 950 }
951 }
952 else if(probTPC[7] > 0.6){ // He3
587d006a 953 x[2] = probTPC[7];
954 x3[2] = x[2];
78a615f6 955 if(TMath::Abs(nsigmaTPC[7]) < 5){
587d006a 956 if(fV2) contV0[iV0]->Fill(7,aodTrack->Pt()*2,v2V0,x);
a72e0394 957 if(fV3) contV0v3[iV0]->Fill(7,aodTrack->Pt()*2,v3V0,x3);
afa8df58 958 }
afa8df58 959 }
960 } // end side loop
961 } // end track loop
962
963 // Fill EP distribution histograms
964 if(fV2) fPhiRPv0A->Fill(iC,evPlAngV0ACor2);
965 if(fV2) fPhiRPv0C->Fill(iC,evPlAngV0CCor2);
966
967 if(fV3) fPhiRPv0Av3->Fill(iC,evPlAngV0ACor3);
968 if(fV3) fPhiRPv0Cv3->Fill(iC,evPlAngV0CCor3);
969
970 // TPC EP needed for resolution studies (TPC subevent)
971 Double_t Qx2 = 0, Qy2 = 0;
972 Double_t Qx3 = 0, Qy3 = 0;
973
974 for(Int_t iT = 0; iT < nAODTracks; iT++) {
975
976 AliAODTrack* aodTrack = aodEvent->GetTrack(iT);
977
978 if (!aodTrack){
979 aodTrack->Delete();
980 continue;
981 }
982
983 Bool_t trkFlag = aodTrack->TestFilterBit(1);
984
985 if ((TMath::Abs(aodTrack->Eta()) > 0.8) || (aodTrack->Pt() < 0.2) || (aodTrack->GetTPCNcls() < 70) || !trkFlag)
986 continue;
987
988 Double_t b[2] = {-99., -99.};
989 Double_t bCov[3] = {-99., -99., -99.};
243fbce7 990 if (!aodTrack->PropagateToDCA(fOutputAOD->GetPrimaryVertex(), fOutputAOD->GetMagneticField(), 100., b, bCov))
afa8df58 991 continue;
992
993 if ((TMath::Abs(b[0]) > 3.0) || (TMath::Abs(b[1]) > 2.4))
994 continue;
995
996 Qx2 += TMath::Cos(2*aodTrack->Phi());
997 Qy2 += TMath::Sin(2*aodTrack->Phi());
998 Qx3 += TMath::Cos(3*aodTrack->Phi());
999 Qy3 += TMath::Sin(3*aodTrack->Phi());
1000
1001 }
1002
1003 evPlAng2 = TMath::ATan2(Qy2, Qx2)/2.;
1004 evPlAng3 = TMath::ATan2(Qy3, Qx3)/3.;
1005
1006 // Fill histograms needed for resolution evaluation
1007 if(fV2) fHResTPCv0A2->Fill(Double_t(iC), TMath::Cos(2*(evPlAng2 - evPlAngV0ACor2)));
1008 if(fV2) fHResTPCv0C2->Fill(Double_t(iC), TMath::Cos(2*(evPlAng2 - evPlAngV0CCor2)));
1009 if(fV2) fHResv0Cv0A2->Fill(Double_t(iC), TMath::Cos(2*(evPlAngV0ACor2 - evPlAngV0CCor2)));
1010
1011 if(fV3) fHResTPCv0A3->Fill(Double_t(iC), TMath::Cos(3*(evPlAng3 - evPlAngV0ACor3)));
1012 if(fV3) fHResTPCv0C3->Fill(Double_t(iC), TMath::Cos(3*(evPlAng3 - evPlAngV0CCor3)));
1013 if(fV3) fHResv0Cv0A3->Fill(Double_t(iC), TMath::Cos(3*(evPlAngV0ACor3 - evPlAngV0CCor3)));
1014 }
1015
1016}
1017
1018//_____________________________________________________________________________
1019Float_t AliAnalysisTaskVnV0::GetVertex(AliAODEvent* aod) const
1020{
1021
1022 Float_t zvtx = -999;
1023
1024 const AliAODVertex* vtxAOD = aod->GetPrimaryVertex();
1025 if (!vtxAOD)
1026 return zvtx;
1027 if(vtxAOD->GetNContributors()>0)
1028 zvtx = vtxAOD->GetZ();
1029
1030 return zvtx;
1031}
1032//_____________________________________________________________________________
1033void AliAnalysisTaskVnV0::Terminate(Option_t *)
1034{
1035 // Terminate loop
1036 Printf("Terminate()");
1037}
1038//_____________________________________________________________________________
1039void AliAnalysisTaskVnV0::OpenInfoCalbration(Int_t run){
1040 TString oadbfilename = "$ALICE_ROOT/OADB/PWGCF/VZERO/VZEROcalibEP.root";
1041 TFile *foadb = TFile::Open(oadbfilename.Data());
1042
1043 if(!foadb){
1044 printf("OADB file %s cannot be opened\n",oadbfilename.Data());
1045 return;
1046 }
1047
1048 AliOADBContainer *cont = (AliOADBContainer*) foadb->Get("hMultV0BefCorr");
1049 if(!cont){
1050 printf("OADB object hMultV0BefCorr is not available in the file\n");
1051 return;
1052 }
1053
1054 if(!(cont->GetObject(run))){
1055 printf("OADB object hMultV0BefCorr is not available for run %i (used run 137366)\n",run);
1056 run = 137366;
1057 }
1058 fMultV0 = ((TH2F *) cont->GetObject(run))->ProfileX();
1059
1060 TF1 *fpol0 = new TF1("fpol0","pol0");
1061 fMultV0->Fit(fpol0,"","",0,31);
1062 fV0Cpol = fpol0->GetParameter(0);
1063 fMultV0->Fit(fpol0,"","",32,64);
1064 fV0Apol = fpol0->GetParameter(0);
1065
1066 for(Int_t iside=0;iside<2;iside++){
1067 for(Int_t icoord=0;icoord<2;icoord++){
1068 for(Int_t i=0;i < nCentrBin;i++){
1069 char namecont[100];
1070 if(iside==0 && icoord==0)
d6a1c304 1071 snprintf(namecont,100,"hQxc2_%i",i);
afa8df58 1072 else if(iside==1 && icoord==0)
d6a1c304 1073 snprintf(namecont,100,"hQxa2_%i",i);
afa8df58 1074 else if(iside==0 && icoord==1)
d6a1c304 1075 snprintf(namecont,100,"hQyc2_%i",i);
afa8df58 1076 else if(iside==1 && icoord==1)
d6a1c304 1077 snprintf(namecont,100,"hQya2_%i",i);
afa8df58 1078
1079 cont = (AliOADBContainer*) foadb->Get(namecont);
1080 if(!cont){
1081 printf("OADB object %s is not available in the file\n",namecont);
1082 return;
1083 }
1084
1085 if(!(cont->GetObject(run))){
1086 printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
1087 run = 137366;
1088 }
1089 fMeanQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
1090 fWidthQ[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
1091
1092 //for v3
1093 if(iside==0 && icoord==0)
d6a1c304 1094 snprintf(namecont,100,"hQxc3_%i",i);
afa8df58 1095 else if(iside==1 && icoord==0)
d6a1c304 1096 snprintf(namecont,100,"hQxa3_%i",i);
afa8df58 1097 else if(iside==0 && icoord==1)
d6a1c304 1098 snprintf(namecont,100,"hQyc3_%i",i);
afa8df58 1099 else if(iside==1 && icoord==1)
d6a1c304 1100 snprintf(namecont,100,"hQya3_%i",i);
afa8df58 1101
1102 cont = (AliOADBContainer*) foadb->Get(namecont);
1103 if(!cont){
1104 printf("OADB object %s is not available in the file\n",namecont);
1105 return;
1106 }
1107
1108 if(!(cont->GetObject(run))){
1109 printf("OADB object %s is not available for run %i (used run 137366)\n",namecont,run);
1110 run = 137366;
1111 }
1112 fMeanQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetMean();
1113 fWidthQv3[i][iside][icoord] = ((TH1F *) cont->GetObject(run))->GetRMS();
1114
1115 }
1116 }
1117 }
1118}