AliPhysics  6f1d526 (6f1d526)
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
AliAnalysisTaskVnZDC.cxx
Go to the documentation of this file.
1 /*************************************************************************
2  * Copyright(c) 1998-2008, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
17 // AliAnalysisTaskVnZDC:
18 // Author: Rihan Haque (mhaque@cern.ch)
20 
21 
22 #include "Riostream.h" //needed as include
23 
24 class TFile;
25 class TList;
26 class AliAnalysisTaskSE;
27 
28 #include "TProfile.h" //needed as include
29 #include "TNtuple.h"
30 #include "TProfile2D.h"
31 #include "TList.h"
32 #include "TH1.h"
33 #include "TH2.h"
34 #include "TH3.h"
35 #include "AliMultSelection.h"
36 #include "AliVVertex.h"
37 
38 #include "AliAnalysisManager.h"
39 #include "AliFlowEventSimple.h"
40 #include "AliAnalysisTaskVnZDC.h"
41 #include "AliFlowCommonHist.h"
43 
44 #include "AliAnalysisManager.h"
45 #include "AliInputEventHandler.h"
46 #include "AliVEvent.h"
47 #include "AliESD.h"
48 #include "AliESDEvent.h"
49 #include "AliESDHeader.h"
50 #include "AliESDInputHandler.h"
51 #include "AliESDZDC.h"
52 #include "AliMultiplicity.h"
53 #include "AliAnalysisUtils.h"
54 #include "AliAODHandler.h"
55 #include "AliAODTrack.h"
56 #include "AliAODEvent.h"
57 #include "AliAODHeader.h"
58 #include "AliAODVertex.h"
59 #include "AliAODVZERO.h"
60 #include "AliAODZDC.h"
61 #include "AliAODMCHeader.h"
62 #include "AliMCEventHandler.h"
63 #include "AliMCEvent.h"
64 #include "AliHeader.h"
65 #include "AliVParticle.h"
66 #include "AliStack.h"
67 #include "AliAODMCParticle.h"
68 #include "AliAnalysisTaskSE.h"
69 #include "AliGenEventHeader.h"
70 #include "AliPhysicsSelectionTask.h"
71 #include "AliPhysicsSelection.h"
72 #include "AliBackgroundSelection.h"
73 #include "AliTriggerAnalysis.h"
74 #include "AliCentrality.h"
75 #include "AliLog.h"
76 
77 using std::endl;
78 using std::cout;
79 
80 
82 
83 //________________________________________________________________________
84 AliAnalysisTaskVnZDC::AliAnalysisTaskVnZDC(const char *name, Bool_t usePhiWeights) :
85 AliAnalysisTaskSE(name),
86 fEvent(NULL),
87 fListHistos(NULL),
88 fMinimalBook(kFALSE),
89 fUsePhiWeights(usePhiWeights),
90 fListWeights(NULL),
91 fRelDiffMsub(1.0),
92 fApplyCorrectionForNUA(kFALSE),
93 fHarmonic(2),
94 fNormalizationType(1),
95 fNCentBins(9),
96 fTotalQvector(NULL),
97 fievent(0),
98 frunflag(0),
99 fHist_Vertex_XY(NULL),
100 fHist_Vx_vs_runnum(NULL),
101 fHist_Vy_vs_runnum(NULL),
102 fHist_Vz_vs_runnum(NULL),
103 fHist_tracks_vs_runnum(NULL),
104 fHist_Vx_ArrayFinder(NULL),
105 fHist_Vy_ArrayFinder(NULL),
106 fHist_Vz_ArrayFinder(NULL),
107 fHist_Vertex_Vz(NULL),
108 fHist_Event_count(NULL),
109 fHist_Cent_count1(NULL),
110 fHist_Cent_count2(NULL),
111 fHist_EventperRun(NULL),
112 fHist_Psi1_zdnA(NULL),
113 fHist_Psi1_zdnC(NULL),
114 fHist_Psi1_zdnA_after1(NULL),
115 fHist_Psi1_zdnC_after1(NULL),
116 fHist_Psi1_zdnA_after2(NULL),
117 fHist_Psi1_zdnC_after2(NULL),
118 fHist_vBincount(NULL),
119 fZDCESEList(NULL),
120 fZListDummy(NULL),
121 fFBEffiList1(NULL),
122 fPileUpMultSelCount(NULL),
123 fPileUpCount(NULL),
124 fMultSelection(NULL),
125 fAnalysisUtil(NULL),
126 fHist_v2xV1_ZDN_Norm_cosXX(NULL),
127 fHist_v2xV1_ZDN_Refm_cosXX(NULL),
128 fHist_v2xV1_ZDN_Cent_cosXX(NULL),
129 fHist_v2xV1_ZDN_Norm_cosYY(NULL),
130 fHist_v2xV1_ZDN_Refm_cosYY(NULL),
131 fHist_v2xV1_ZDN_Cent_cosYY(NULL),
132 fHist_v2xV1_ZDN_Norm_sinXY(NULL),
133 fHist_v2xV1_ZDN_Refm_sinXY(NULL),
134 fHist_v2xV1_ZDN_Cent_sinXY(NULL),
135 fHist_v2xV1_ZDN_Norm_sinYX(NULL),
136 fHist_v2xV1_ZDN_Refm_sinYX(NULL),
137 fHist_v2xV1_ZDN_Cent_sinYX(NULL),
138 fHist_v2xV1_ZDN_Norm_All(NULL),
139 fHist_v2xV1_ZDN_Refm_All(NULL),
140 fHist_v2xV1_ZDN_Cent_All(NULL),
141 fHist_ZDN_resol_Norm_All(NULL),
142 fHist_ZDN_resol_Cent_All(NULL),
143 fHist_ZDN_resol_Refm_All(NULL),
144 fHist_ZDN_resol_Norm_cos(NULL),
145 fHist_ZDN_resol_Refm_cos(NULL),
146 fHist_ZDN_resol_Cent_cos(NULL),
147 fHist_ZDN_resol_Norm_sin(NULL),
148 fHist_ZDN_resol_Refm_sin(NULL),
149 fHist_ZDN_resol_Cent_sin(NULL),
150 fHist_ZDCn_A_XYvsRun(NULL),
151 fHist_ZDCn_C_XYvsRun(NULL),
152 fAvg_Cent_vs_Vz_Cent_woCut(NULL),
153 fAvg_Cent_vs_Vz_Peri_woCut(NULL),
154 fAvg_Cent_vs_Vx_Cent_woCut(NULL),
155 fAvg_Cent_vs_Vx_Peri_woCut(NULL),
156 fAvg_Cent_vs_Vy_Cent_woCut(NULL),
157 fAvg_Cent_vs_Vy_Peri_woCut(NULL),
158 fAvg_Cent_vs_Vz_Cent_wCuts(NULL),
159 fAvg_Cent_vs_Vz_Peri_wCuts(NULL),
160 fAvg_Cent_vs_Vx_Cent_wCuts(NULL),
161 fAvg_Cent_vs_Vx_Peri_wCuts(NULL),
162 fAvg_Cent_vs_Vy_Cent_wCuts(NULL),
163 fAvg_Cent_vs_Vy_Peri_wCuts(NULL),
164 fTPCV0M_CentDiff_vs_Vz(NULL),
165 fTPCV0M_CentDiff_vs_Vx(NULL),
166 fTPCV0M_CentDiff_vs_Vy(NULL),
167 fWeight_Cent(NULL),
168 fCent_fromDATA(NULL),
169 fRejectPileUpTight(kFALSE),
170 fRejectPileUp(kFALSE),
171 checkOnce(0),
172 vxBin(0),
173 vyBin(0),
174 vzBin(0)
175 {
176  for(int i=0;i<90;i++)
177  runNums[i] = 0;
178 
179  for(int i=0;i<90;i++){
180  for(int j=0;j<10;j++){
181  fHist_znCx_V0_VxVy[i][j] = NULL;
182  fHist_znCy_V0_VxVy[i][j] = NULL;
183  fHist_znAx_V0_VxVy[i][j] = NULL;
184  fHist_znAy_V0_VxVy[i][j] = NULL;
185  }
186  fHist_ZDCA_En_Run[i] = NULL;
187  fHist_ZDCC_En_Run[i] = NULL;
188  }
189 
190  for(int i=0;i<2;i++)
191  {
192  VxCut[i] = 0;
193  VyCut[i] = 0;
194  VzCut[i] = 0;
195  }
196 
197  for(int i=0;i<4;i++){
198  for(int j=0;j<5;j++){
199  fHist_X_vs_Obs_before[i][j] = NULL;
200  fHist_X_vs_Obs_after1[i][j] = NULL;
201  fHist_XX_vs_Obs_before[i][j] = NULL;
202  fHist_XX_vs_Obs_after1[i][j] = NULL;
203  }
204  }
205 
206  for(int i=0;i<10;i++){
207  fFB_Efficiency_Cent[i] = NULL;
208  }
209 
210 
211  DefineInput(1, AliFlowEventSimple::Class()); // Input slot #0 works with an AliFlowEventSimple
212  DefineOutput(1,TList::Class());
213  DefineOutput(2,TList::Class());
214 
215  fDataSet="2010";
216  fAnalysisSet="recenter1";
217 
218  fTotalQvector = new TString("QaQb"); // Total Q-vector is: "QaQb" (means Qa+Qb), "Qa" or "Qb"
219  //AliDebug(2,
220  cout<<"AliAnalysisTaskVnZDC::Constructor is called..!! fAnalysisSet = "<<fAnalysisSet.Data()<<endl;
221 
222 }//-------------constructor-----------------
223 
224 //________________________________________________
227 fEvent(NULL),
228 fListHistos(NULL),
229 fMinimalBook(kFALSE),
230 fUsePhiWeights(kFALSE),
231 fListWeights(NULL),
232 fRelDiffMsub(1.0),
233 fApplyCorrectionForNUA(kFALSE),
234 fHarmonic(2),
235 fNormalizationType(1),
236 fNCentBins(9),
237 fTotalQvector(NULL),
238 fievent(0),
239 frunflag(0),
240 fHist_Vertex_XY(NULL),
241 fHist_Vx_vs_runnum(NULL),
242 fHist_Vy_vs_runnum(NULL),
243 fHist_Vz_vs_runnum(NULL),
244 fHist_tracks_vs_runnum(NULL),
245 fHist_Vx_ArrayFinder(NULL),
246 fHist_Vy_ArrayFinder(NULL),
247 fHist_Vz_ArrayFinder(NULL),
248 fHist_Vertex_Vz(NULL),
249 fHist_Event_count(NULL),
250 fHist_Cent_count1(NULL),
251 fHist_Cent_count2(NULL),
252 fHist_EventperRun(NULL),
253 fHist_Psi1_zdnA(NULL),
254 fHist_Psi1_zdnC(NULL),
255 fHist_Psi1_zdnA_after1(NULL),
256 fHist_Psi1_zdnC_after1(NULL),
257 fHist_Psi1_zdnA_after2(NULL),
258 fHist_Psi1_zdnC_after2(NULL),
259 fHist_vBincount(NULL),
260 fZDCESEList(NULL),
261 fZListDummy(NULL),
262 fFBEffiList1(NULL),
263 fPileUpMultSelCount(NULL),
264 fPileUpCount(NULL),
265 fMultSelection(NULL),
266 fAnalysisUtil(NULL),
267 fHist_v2xV1_ZDN_Norm_cosXX(NULL),
268 fHist_v2xV1_ZDN_Refm_cosXX(NULL),
269 fHist_v2xV1_ZDN_Cent_cosXX(NULL),
270 fHist_v2xV1_ZDN_Norm_cosYY(NULL),
271 fHist_v2xV1_ZDN_Refm_cosYY(NULL),
272 fHist_v2xV1_ZDN_Cent_cosYY(NULL),
273 fHist_v2xV1_ZDN_Norm_sinXY(NULL),
274 fHist_v2xV1_ZDN_Refm_sinXY(NULL),
275 fHist_v2xV1_ZDN_Cent_sinXY(NULL),
276 fHist_v2xV1_ZDN_Norm_sinYX(NULL),
277 fHist_v2xV1_ZDN_Refm_sinYX(NULL),
278 fHist_v2xV1_ZDN_Cent_sinYX(NULL),
279 fHist_ZDN_resol_Norm_cos(NULL),
280 fHist_ZDN_resol_Refm_cos(NULL),
281 fHist_ZDN_resol_Cent_cos(NULL),
282 fHist_ZDN_resol_Norm_sin(NULL),
283 fHist_ZDN_resol_Refm_sin(NULL),
284 fHist_ZDN_resol_Cent_sin(NULL),
285 fHist_v2xV1_ZDN_Norm_All(NULL),
286 fHist_v2xV1_ZDN_Refm_All(NULL),
287 fHist_v2xV1_ZDN_Cent_All(NULL),
288 fHist_ZDN_resol_Norm_All(NULL),
289 fHist_ZDN_resol_Cent_All(NULL),
290 fHist_ZDN_resol_Refm_All(NULL),
291 fHist_ZDCn_A_XYvsRun(NULL),
292 fHist_ZDCn_C_XYvsRun(NULL),
293 fAvg_Cent_vs_Vz_Cent_woCut(NULL),
294 fAvg_Cent_vs_Vz_Peri_woCut(NULL),
295 fAvg_Cent_vs_Vx_Cent_woCut(NULL),
296 fAvg_Cent_vs_Vx_Peri_woCut(NULL),
297 fAvg_Cent_vs_Vy_Cent_woCut(NULL),
298 fAvg_Cent_vs_Vy_Peri_woCut(NULL),
299 fAvg_Cent_vs_Vz_Cent_wCuts(NULL),
300 fAvg_Cent_vs_Vz_Peri_wCuts(NULL),
301 fAvg_Cent_vs_Vx_Cent_wCuts(NULL),
302 fAvg_Cent_vs_Vx_Peri_wCuts(NULL),
303 fAvg_Cent_vs_Vy_Cent_wCuts(NULL),
304 fAvg_Cent_vs_Vy_Peri_wCuts(NULL),
305 fTPCV0M_CentDiff_vs_Vz(NULL),
306 fTPCV0M_CentDiff_vs_Vx(NULL),
307 fTPCV0M_CentDiff_vs_Vy(NULL),
308 fWeight_Cent(NULL),
309 fCent_fromDATA(NULL),
310 fRejectPileUpTight(kFALSE),
311 fRejectPileUp(kFALSE),
312 checkOnce(0),
313 vxBin(0),
314 vyBin(0),
315 vzBin(0)
316 {
317  for(int i=0;i<90;i++)
318  runNums[i] = 0;
319 
320  for(int i=0;i<90;i++){
321  for(int j=0;j<10;j++){
322  fHist_znCx_V0_VxVy[i][j] = NULL;
323  fHist_znCy_V0_VxVy[i][j] = NULL;
324  fHist_znAx_V0_VxVy[i][j] = NULL;
325  fHist_znAy_V0_VxVy[i][j] = NULL;
326  }
327  fHist_ZDCA_En_Run[i] = NULL;
328  fHist_ZDCC_En_Run[i] = NULL;
329  }
330 
331  for(int i=0;i<2;i++)
332  {
333  VxCut[i] = 0;
334  VyCut[i] = 0;
335  VzCut[i] = 0;
336  }
337 
338  for(int i=0;i<4;i++){
339  for(int j=0;j<5;j++){
340  fHist_X_vs_Obs_before[i][j] = NULL;
341  fHist_X_vs_Obs_after1[i][j] = NULL;
342  fHist_XX_vs_Obs_before[i][j] = NULL;
343  fHist_XX_vs_Obs_after1[i][j] = NULL;
344  }
345  }
346 
347  for(int i=0;i<10;i++){
348  fFB_Efficiency_Cent[i] = NULL;
349  }
350 
351 
352  fDataSet="2010";
353  fAnalysisSet="recenter1";
354 
355  // AliDebug(2,"AliAnalysisTaskVnZDC::AliAnalysisTaskVnZDC()");
356 }
357 
358 
359 
360 
361 
362 //________________________________________________________________________
364 {
365  //delete [] fSP; //is this not deleting all fSP[i] ?
366  delete fListHistos;
367  delete fHist_Vx_ArrayFinder;
368  delete fHist_Vy_ArrayFinder;
369  delete fHist_Vz_ArrayFinder;
370  delete fMultSelection;
371  delete fAnalysisUtil;
372  delete fEvent;
373  delete fCent_fromDATA;
374  delete fWeight_Cent;
375  //delete fHist_ZDCn_A_XYvsRun;
376  //delete fHist_ZDCn_C_XYvsRun;
377 
378  if(fAnalysisSet=="recenter2"){
379  for(int i=0;i<90;i++){
380  for(int j=0;j<10;j++){
381  delete fHist_znCx_V0_VxVy[i][j];
382  delete fHist_znCy_V0_VxVy[i][j];
383  delete fHist_znAx_V0_VxVy[i][j];
384  delete fHist_znAy_V0_VxVy[i][j];
385  }
386  delete fHist_ZDCA_En_Run[i];
387  delete fHist_ZDCC_En_Run[i];
388  }
389  }
390 
391 
392  for(int i=0;i<10;i++){
393  delete fFB_Efficiency_Cent[i];
394  }
395 
396  delete fListHistos;
397  delete fFBEffiList1;
398  delete fListWeights;
399  delete fZDCESEList;
400  delete fZListDummy;
401 
402  printf("\n\n ~AliAnalysisTaskVnZDC::Destructor is called..!!\n\n");
403 
404 }
405 
406 //________________________________________________________________________
408 {
409  //printf("\n\n... AliAnalysisTaskVnZDC::CreateOutputObjects() beginning 0 ...\n");
410  //AliDebug(2,"AliAnalysisTaskVnZDC::CreateOutputObjects() is called....");
411 
412  Int_t runArray_2010[89] = {139510, 139507, 139505, 139503, 139465, 139438, 139437, 139360, 139329, 139328, 139314, 139310, 139309, 139173, 139107, 139105, 139038, 139037, 139036, 139029, 139028, 138872, 138871, 138870, 138837, 138732, 138730, 138666, 138662, 138653, 138652, 138638, 138624, 138621, 138583, 138582, 138579, 138578, 138534, 138469, 138442, 138439, 138438, 138396, 138364, 138275, 138225, 138201, 138197, 138192, 138190, 137848, 137844, 137752, 137751, 137724, 137722, 137718, 137704, 137693, 137692, 137691, 137686, 137685, 137639, 137638, 137608, 137595, 137549, 137546, 137544, 137541, 137539, 137531, 137530, 137443, 137441, 137440, 137439, 137434, 137432, 137431, 137243, 137236, 137235, 137232, 137231, 137162, 137161};
413 
414  Int_t runArray_2011[68] = {167915, 168115, 168460, 169035, 169238, 169859, 170228, 167920, 168310, 168464, 169091, 169411, 169923, 170230, 167985, 168311, 168467, 169094, 169415, 170027, 170268, 167987, 168322, 168511, 169138, 169417, 170081, 170269, 167988, 168325, 168512, 169144, 169835, 170155, 170270, 168069, 168341, 168514, 169145, 169837, 170159, 170306, 168076, 168342, 168777, 169148, 169838, 170163, 170308, 168105, 168361, 168826, 169156, 169846, 170193, 170309, 168107, 168362, 168988, 169160, 169855, 170203, 168108, 168458, 168992, 169167, 169858, 170204};
415 
416  Int_t runArray_2015[90] = {246994, 246991, 246989, 246984, 246982, 246980, 246948, 246945, 246928, 246871, 246870, 246867, 246865, 246864, 246859, 246858, 246851, 246847, 246846, 246845, 246844, 246810, 246809, 246808, 246807, 246805, 246804, 246766, 246765, 246763, 246760, 246759, 246758, 246757, 246751, 246750, 246676, 246675, 246540, 246495, 246493, 246488, 246487, 246434, 246431, 246428, 246424, 246276, 246275, 246272, 246271, 246225, 246222, 246217, 246185, 246182, 246181, 246180, 246178, 246153, 246152, 246151, 246148, 246115, 246113, 246089, 246087, 246053, 246052, 246049, 246048, 246042, 246037, 246036, 246012, 246003, 246001, 245963, 245954, 245952, 245949, 245923, 245833, 245831, 245829, 245705, 245702, 245700, 245692, 245683};
417 
418  //Int_t runArray_2015[77] = {246994, 246991, 246989, 246984, 246982, 246948, 246945, 246928, 246851, 246847, 246846, 246845, 246844, 246810, 246809, 246808, 246807, 246805, 246804, 246766, 246765, 246763, 246760, 246759, 246758, 246757, 246751, 246750, 246495, 246493, 246488, 246487, 246434, 246431, 246424, 246276, 246275, 246272, 246271, 246225, 246222, 246217, 246185, 246182, 246181, 246180, 246178, 246153, 246152, 246151, 246148, 246115, 246113, 246089, 246087, 246053, 246052, 246049, 246048, 246042, 246037, 246036, 246012, 246003, 246001, 245963, 245954, 245952, 245949, 245923, 245833, 245831, 245829, 245705, 245702, 245692, 245683};
419 
420 
421  if(fDataSet=="2010"){
422  frunflag = 89;
423  for(int i=0;i<frunflag;i++)
424  runNums[i] = runArray_2010[i];
425  }
426  if(fDataSet=="2011"){
427  frunflag = 10; //<--------- 2011 runnumbers have to be checked and put...
428  for(int i=0;i<frunflag;i++)
429  runNums[i] = runArray_2011[i];
430  }
431  if(fDataSet=="2015"){
432  frunflag = 90;
433  for(int i=0;i<frunflag;i++)
434  runNums[i] = runArray_2015[i];
435  }
436 
437 
438  Int_t totalQ = 0;
439  if(fTotalQvector->Contains("Qa")) totalQ += 1;
440  if(fTotalQvector->Contains("Qb")) totalQ += 2;
441 
442  fListHistos = new TList();
443  fListHistos->SetOwner(kTRUE);
444 
445 
446  fHist_Event_count = new TH1F("fHist_Event_count"," ",20,0,20);
448  fHist_vBincount = new TH1F("fHist_vBincount"," ",500,0,500);
450  fHist_EventperRun = new TH1F("fHist_EventperRun","", frunflag,0,frunflag);
452  fHist_Psi1_zdnA = new TH1F("fHist_Psi1_zdnA","", 200, 0,2.*TMath::Pi());
454  fHist_Psi1_zdnC = new TH1F("fHist_Psi1_zdnC","", 200, 0,2.*TMath::Pi());
456  fHist_Psi1_zdnA_after1 = new TH1F("fHist_Psi1_zdnA_after1","", 200, 0,2.*TMath::Pi());
458  fHist_Psi1_zdnC_after1 = new TH1F("fHist_Psi1_zdnC_after1","", 200, 0,2.*TMath::Pi());
460  fHist_Psi1_zdnA_after2 = new TH1F("fHist_Psi1_zdnA_after2","", 200, 0,2.*TMath::Pi());
462  fHist_Psi1_zdnC_after2 = new TH1F("fHist_Psi1_zdnC_after2","", 200, 0,2.*TMath::Pi());
464 
465 
466  fHist_Vertex_Vz = new TH1F("fHist_Vertex_Vz_dist","Vertex Z Distribution",200, -15, 15);
468  fHist_Vertex_XY = new TH2F("fHist_Vertex_XY_dist","Vertex XY Distributn.",100,-0.2,0.2,100,0.0,0.4);
470 
471  // Vx,Vy,Vz vs runnumber:
472  fHist_Vx_vs_runnum = new TProfile("fHist_Vx_vs_runnum","<Vx>_vs_runnum",frunflag,0,frunflag);
474  fHist_Vy_vs_runnum = new TProfile("fHist_Vy_vs_runnum","<Vy>_vs_runnum",frunflag,0,frunflag);
476  fHist_Vz_vs_runnum = new TProfile("fHist_Vz_vs_runnum","<Vz>_vs_runnum",frunflag,0,frunflag);
478  fHist_tracks_vs_runnum = new TProfile("fHist_tracks_vs_runnum","<nTracks>_vs_runnum",frunflag,0,frunflag);
480 
481  fHist_ZDCn_A_XYvsRun = new TH3F("fHist_ZDCn_A_XYvsRun","ZDC XY vs run",80,-2.0,2.0,80,-2.0,2.0,89,0,89);
483  fHist_ZDCn_C_XYvsRun = new TH3F("fHist_ZDCn_C_XYvsRun","ZDC XY vs run",80,-2.0,2.0,80,-2.0,2.0,89,0,89);
485 
486 
487  Double_t centRange[12] = {0,5,10,20,30,40,50,60,70,80,90,100};
488 
489  fPileUpCount = new TH1F("fPileUpCount", "fPileUpCount", 9, 0., 9.);
490  fPileUpCount->GetXaxis()->SetBinLabel(1,"plpMV");
491  fPileUpCount->GetXaxis()->SetBinLabel(2,"fromSPD");
492  fPileUpCount->GetXaxis()->SetBinLabel(3,"RefMultiplicityComb08");
493  fPileUpCount->GetXaxis()->SetBinLabel(4,"IncompleteDAQ");
494  fPileUpCount->GetXaxis()->SetBinLabel(5,"abs(V0M-CL1)>7.5");
495  fPileUpCount->GetXaxis()->SetBinLabel(6,"missingVtx");
496  fPileUpCount->GetXaxis()->SetBinLabel(7,"inconsistentVtx");
497  fPileUpCount->GetXaxis()->SetBinLabel(8,"multESDTPCDif");
498  fPileUpCount->GetXaxis()->SetBinLabel(9,"extraPileUpMultSel");
500 
501  fPileUpMultSelCount = new TH1F("fPileUpMultSelCount", "fPileUpMultSelCount", 8, 0., 8.);
502  fPileUpMultSelCount->GetXaxis()->SetBinLabel(1,"IsNotPileup");
503  fPileUpMultSelCount->GetXaxis()->SetBinLabel(2,"IsNotPileupMV");
504  fPileUpMultSelCount->GetXaxis()->SetBinLabel(3,"IsNotPileupInMultBins");
505  fPileUpMultSelCount->GetXaxis()->SetBinLabel(4,"InconsistentVertices");
506  fPileUpMultSelCount->GetXaxis()->SetBinLabel(5,"TrackletVsCluster");
507  fPileUpMultSelCount->GetXaxis()->SetBinLabel(6,"AsymmetricInVZERO");
508  fPileUpMultSelCount->GetXaxis()->SetBinLabel(7,"IncompleteDAQ");
509  fPileUpMultSelCount->GetXaxis()->SetBinLabel(8,"GoodVertex2016");
511 
512 
513 
514 
515  //-------- define and add the recenter histograms ----------------
516 if(fDataSet=="2010") {
517  vxBin = 8;
518  vyBin = 8;
519  //vxBin = 20;
520  //vyBin = 20;
521  vzBin = 10; // j<10
522 
523  VxCut[0] = -0.035;
524  VxCut[1] = 0.020;
525  VyCut[0] = 0.144;
526  VyCut[1] = 0.214;
527  VzCut[0] = -10.0;
528  VzCut[1] = 10.0;
529  }
530 
531  if(fDataSet=="2015"){
532  vxBin = 8;
533  vyBin = 8;
534  //vxBin = 20;
535  //vyBin = 20;
536 
537  vzBin = 10; // j<10
538 
539  VxCut[0] = 0.060;
540  VxCut[1] = 0.086;
541  VyCut[0] = 0.321;
542  VyCut[1] = 0.345;
543  VzCut[0] = -10.0;
544  VzCut[1] = 10.0;
545  }
546  if(fDataSet=="2011"){
547  AliDebug(2,"\n\n !!** WARNING ***!! \n cuts not defined for DATASET: %s\n ...EXIT...\n\n)");
548  exit(0);
549  }
550 
551 
552  const int NbinVt = (vyBin-1)*vxBin + vxBin;
553 
554  fHist_Vx_ArrayFinder = new TH1F("fHist_Vx_ArrayFinder","",vxBin,VxCut[0],VxCut[1]);
555  fHist_Vy_ArrayFinder = new TH1F("fHist_Vy_ArrayFinder","",vyBin,VyCut[0],VyCut[1]);
556  fHist_Vz_ArrayFinder = new TH1F("fHist_Vz_ArrayFinder","",vzBin,VzCut[0],VzCut[1]);
557 
558 
559 
560 
561 
562  if(fAnalysisSet=="recenter2") {
563 
564  TString sNameQn[4] = {"Xa","Xc","Ya","Yc"};
565  TString sNameQn2[4] = {"XaXc","YaYc","XcYa","YcXa"};
566  TString sNameVs[5] = {"Cent","Mult","Vx","Vy","Vz"}; // ********** to add vs cent ****** !!!!!!!!!!!!
567  Int_t nBinNumVs[5] = {100, 500, 40, 40,400}; // number of bins for "Mult","Vx","Vy","Vz"
568  Double_t lBinLowVs[5] = { 0, 0,-0.1,0.1,-10}; // bin low value: "Mult","Vx","Vy","Vz"
569  Double_t lBinHighVs[5] = {100,2000, 0.1,0.3, 10};// bin high value: "Mult","Vx","Vy","Vz"
570 
571  char name[100];
572 
573 
574  for(int i=0;i<4;i++) {
575  for(int j=0;j<5;j++) {
576  //store: X,Y position for recenter:
577  sprintf(name,"fHist_%s_vs_%s_before",static_cast<const char*>(sNameQn[i]),static_cast<const char*>(sNameVs[j]));
578  fHist_X_vs_Obs_before[i][j] = new TProfile(name,"", nBinNumVs[j], lBinLowVs[j], lBinHighVs[j],"");
580 
581  sprintf(name,"fHist_%s_vs_%s_after1",static_cast<const char*>(sNameQn[i]),static_cast<const char*>(sNameVs[j]));
582  fHist_X_vs_Obs_after1[i][j] = new TProfile(name,"", nBinNumVs[j], lBinLowVs[j], lBinHighVs[j],"");
584 
585  sprintf(name,"fHist_%s_vs_%s_before",static_cast<const char*>(sNameQn2[i]),static_cast<const char*>(sNameVs[j]));
586  fHist_XX_vs_Obs_before[i][j] = new TProfile(name,"", nBinNumVs[j], lBinLowVs[j], lBinHighVs[j],"");
588 
589  sprintf(name,"fHist_%s_vs_%s_after1",static_cast<const char*>(sNameQn2[i]),static_cast<const char*>(sNameVs[j]));
590  fHist_XX_vs_Obs_after1[i][j] = new TProfile(name,"", nBinNumVs[j], lBinLowVs[j], lBinHighVs[j],"");
592  }
593  }
594 
595  Double_t centRange[12] = {0,5,10,20,30,40,50,60,70,80,90,100};
596 
597  fHist_v2xV1_ZDN_Norm_cosXX = new TProfile("fHist_v2xV1_ZDN_Norm_cosXX","v2 X V1^2 (ZDC) cos(2p)XX",11,centRange);
600  fHist_v2xV1_ZDN_Refm_cosXX = new TProfile("fHist_v2xV1_ZDN_Refm_cosXX","v2 X V1^2 (ZDC) cos(2p)XX", 40,0, 4000);
603  fHist_v2xV1_ZDN_Cent_cosXX = new TProfile("fHist_v2xV1_ZDN_Cent_cosXX","v2 X V1^2 (ZDC) cos(2p)XX",100,0, 100);
606 
607  fHist_v2xV1_ZDN_Norm_cosYY = new TProfile("fHist_v2xV1_ZDN_Norm_cosYY","v2 X V1^2 (ZDC) cos(2p)YY",11,centRange);
610  fHist_v2xV1_ZDN_Refm_cosYY = new TProfile("fHist_v2xV1_ZDN_Refm_cosYY","v2 X V1^2 (ZDC) cos(2p)YY", 40,0, 4000);
613  fHist_v2xV1_ZDN_Cent_cosYY = new TProfile("fHist_v2xV1_ZDN_Cent_cosYY","v2 X V1^2 (ZDC) cos(2p)YY",100,0, 100);
616 
617  fHist_v2xV1_ZDN_Norm_sinXY = new TProfile("fHist_v2xV1_ZDN_Norm_sinXY","v2 X V1^2 (ZDC) sin(2p)XY",11,centRange);
620  fHist_v2xV1_ZDN_Refm_sinXY = new TProfile("fHist_v2xV1_ZDN_Refm_sinXY"," v2 X V1^2 (ZDC) sin(2p)XY", 40,0, 4000);
623  fHist_v2xV1_ZDN_Cent_sinXY = new TProfile("fHist_v2xV1_ZDN_Cent_sinXY"," v2 X V1^2 (ZDC) sin(2p)XY",100,0, 100);
626 
627  fHist_v2xV1_ZDN_Norm_sinYX = new TProfile("fHist_v2xV1_ZDN_Norm_sinYX","v2 X V1^2 (ZDC) sin(2p)YX",11,centRange);
630  fHist_v2xV1_ZDN_Refm_sinYX = new TProfile("fHist_v2xV1_ZDN_Refm_sinYX"," v2 X V1^2 (ZDC) sin(2p)YX", 40,0, 4000);
633  fHist_v2xV1_ZDN_Cent_sinYX = new TProfile("fHist_v2xV1_ZDN_Cent_sinYX"," v2 X V1^2 (ZDC) sin(2p)YX",100,0, 100);
636 
637  fHist_ZDN_resol_Norm_cos = new TProfile("fHist_ZDN_resol_Norm_XX","Resol. V1^2(ZDC) XX",11,centRange);
638  fHist_ZDN_resol_Norm_cos->Sumw2();
640  fHist_ZDN_resol_Refm_cos = new TProfile("fHist_ZDN_resol_Refm_XX","Resol. V1^2(ZDC) XX", 40,0, 4000);
641  fHist_ZDN_resol_Refm_cos->Sumw2();
643  fHist_ZDN_resol_Cent_cos = new TProfile("fHist_ZDN_resol_Cent_XX","Resol. V1^2(ZDC) XX",100,0, 100);
644  fHist_ZDN_resol_Cent_cos->Sumw2();
646 
647  fHist_ZDN_resol_Norm_sin = new TProfile("fHist_ZDN_resol_Norm_YY","Resol. V1^2(ZDC) YY",11,centRange);
648  fHist_ZDN_resol_Norm_sin->Sumw2();
650  fHist_ZDN_resol_Refm_sin = new TProfile("fHist_ZDN_resol_Refm_YY","Resol. V1^2(ZDC) YY",40,0, 4000);
651  fHist_ZDN_resol_Refm_sin->Sumw2();
653  fHist_ZDN_resol_Cent_sin = new TProfile("fHist_ZDN_resol_Cent_YY","Resol. V1^2(ZDC) YY",100,0, 100);
654  fHist_ZDN_resol_Cent_sin->Sumw2();
656 
657  //all terms:
658 
659  fHist_v2xV1_ZDN_Norm_All = new TProfile("fHist_v2xV1_ZDN_Norm_All","v2 X V1^2 (ZDC) all terms",11,centRange);
660  fHist_v2xV1_ZDN_Norm_All->Sumw2();
662  fHist_v2xV1_ZDN_Refm_All = new TProfile("fHist_v2xV1_ZDN_Refm_All","v2 X V1^2 (ZDC) all terms", 40,0, 4000);
663  fHist_v2xV1_ZDN_Refm_All->Sumw2();
665  fHist_v2xV1_ZDN_Cent_All = new TProfile("fHist_v2xV1_ZDN_Cent_All","v2 X V1^2 (ZDC) all terms",100,0, 100);
666  fHist_v2xV1_ZDN_Cent_All->Sumw2();
668 
669  fHist_ZDN_resol_Norm_All = new TProfile("fHist_ZDN_resol_Norm_All","Resol. V1^2(ZDC) All",11,centRange);
670  fHist_ZDN_resol_Norm_All->Sumw2();
672  fHist_ZDN_resol_Refm_All = new TProfile("fHist_ZDN_resol_Refm_All","Resol. V1^2(ZDC) All",40,0, 4000);
673  fHist_ZDN_resol_Refm_All->Sumw2();
675  fHist_ZDN_resol_Cent_All = new TProfile("fHist_ZDN_resol_Cent_All","Resol. V1^2(ZDC) All",100,0, 100);
676  fHist_ZDN_resol_Cent_All->Sumw2();
678 
679 
680  //------------ calculate centrality weight: ------------
681  if(fZDCESEList) {
682  fCent_fromDATA = (TH1F *) fZDCESEList->FindObject("fHist_Cent_afterr_ZDCcut");
683  }
684  else{
685  AliDebug(2,"\n\n ******** Running without Centrality weight !!******** \n\n");
686  fCent_fromDATA = new TH1F("fCent_fromDATA","Centrality distribution",100,0,100);
687  for(int i=1;i<=100;i++){
688  fCent_fromDATA->SetBinContent(i,100);
689  }
690  }
691 
692  Double_t maxCount = -1;
693  Int_t maxCent = -1;
694  Double_t Content,weight,error;
695 
696  for(int i=1;i<=fCent_fromDATA->GetNbinsX();i++){
697  Content = fCent_fromDATA->GetBinContent(i);
698  if(maxCount < Content){
699  maxCount = Content;
700  maxCent = i;
701  }
702  }
703 
704  fWeight_Cent = new TH1F("fWeight_Cent","Weight for centrality",100,0,100);
705 
706  for(int i=1;i<=fCent_fromDATA->GetNbinsX();i++)
707  {
708  Content = fCent_fromDATA->GetBinContent(i);
709  error = fCent_fromDATA->GetBinError(i);
710  weight = maxCount/Content;
711  if(weight>1e9) continue;
712  fWeight_Cent->SetBinContent(i,weight);
713  fWeight_Cent->SetBinError(i,weight/sqrt(Content));
714  //cout<<"cent "<<i-1<<"-"<<i<<" \t wgt = "<<fWeight_Cent->GetBinContent(i)<<" error = "<<fWeight_Cent->GetBinError(i)<<endl;
715  }
716 
718  //------------------------------------------------------
719 
720 
721 
722  //---------Filter bit efficiency----------
723  if(!fFBEffiList1) {
724  printf("\n\n !!** Warning ***!! \n TList FilterBit efficiency not found !!\n\n");
725  }
726 
727  TString name2;
728 
729  if(fFBEffiList1) {
730  for(int i=0;i<10;i++) {
731  fFB_Efficiency_Cent[i] = (TH1D *) fFBEffiList1->FindObject(Form("eff_unbiased_%d",i));
732  name2 = fFB_Efficiency_Cent[i]->GetName();
733  //cout<<i<<" name = "<<name2<<endl;
734  }
735  }
736  else{ // if MC efficiency not found then use weight = 1. Define histograms so that code doesn't crash
737  for(int i=0;i<10;i++){
738  fFB_Efficiency_Cent[i] = new TH1D(Form("eff_unbiased_%d",i),"",100,0,20);
739  for(int j=1;j<=100;j++){
740  fFB_Efficiency_Cent[i]->SetBinContent(j,1.000);
741  }
742  }
743  }
744  //---------Filter bit efficiency----------
745 
746 
747  if(!fZDCESEList) {
748  printf("\n\n !!** ERROR ***!! \n \n TList fZDCESEList not found !!\n .......EXIT...... \n\n)");
749  exit(1);
750  }
751  if(fZDCESEList) {
752  for(int i=0;i<frunflag;i++) {
753  for(int j=0;j<10;j++){
754  fHist_znCx_V0_VxVy[i][j] = (TProfile2D *) fZDCESEList->FindObject(Form("fHist_znCx_V0_Run%d_Vz%d",runNums[i],j+1));
755  fHist_znCy_V0_VxVy[i][j] = (TProfile2D *) fZDCESEList->FindObject(Form("fHist_znCy_V0_Run%d_Vz%d",runNums[i],j+1));
756  fHist_znAx_V0_VxVy[i][j] = (TProfile2D *) fZDCESEList->FindObject(Form("fHist_znAx_V0_Run%d_Vz%d",runNums[i],j+1));
757  fHist_znAy_V0_VxVy[i][j] = (TProfile2D *) fZDCESEList->FindObject(Form("fHist_znAy_V0_Run%d_Vz%d",runNums[i],j+1));
758  if(!fHist_znCx_V0_VxVy[i][j] || !fHist_znCy_V0_VxVy[i][j] || !fHist_znAx_V0_VxVy[i][j] || !fHist_znAy_V0_VxVy[i][j]) {
759  printf("\n\n !!** WARNING *** One/more Recenter1 histograms NOT FOUND\n ...........!! \n\n");
760  //exit(0);
761  }
762  }
763  }
764  }
765  else{ // if recenter file not found then define empty profile histograms so that code doesn't crash
766  //AliDebug(2,"\n\n ******** Running without Recentering 1 !!******** \n\n");
767  printf("\n\n ******** Running without Recentering 1 !!******** \n\n");
768  for(int i=0;i<frunflag;i++){
769  for(int j=0;j<10;j++){
770  fHist_znCx_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znCx_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
771  fHist_znCy_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znCy_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
772  fHist_znAx_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znAx_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
773  fHist_znAy_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znAy_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
774  }
775  }
776  }
777  }//--------------"recenter2"------------------
778 
779 
780  fHist_Cent_count1 = new TH1F("fHist_Cent_before_ZDCcut"," ",100,0,100);
781  fHist_Cent_count2 = new TH1F("fHist_Cent_afterr_ZDCcut"," ",100,0,100);
782 
783 
784  if(fAnalysisSet=="recenter1") {
785 
786  fZDCESEList = new TList();
787  fZDCESEList->SetOwner(kTRUE);
788 
789  fAvg_Cent_vs_Vz_Cent_woCut = new TProfile("fAvg_Cent_vs_Vz_Cent_woCut","<cent> vs Vz",200,-10,10);
791  fAvg_Cent_vs_Vz_Peri_woCut = new TProfile("fAvg_Cent_vs_Vz_Peri_woCut","<peri> vs Vz",200,-10,10);
793  fAvg_Cent_vs_Vx_Cent_woCut = new TProfile("fAvg_Cent_vs_Vx_Cent_woCut","<cent> vs Vx",100,VxCut[0], VxCut[1]);
795  fAvg_Cent_vs_Vx_Peri_woCut = new TProfile("fAvg_Cent_vs_Vx_Peri_woCut","<peri> vs Vx",100,VxCut[0], VxCut[1]);
797  fAvg_Cent_vs_Vy_Cent_woCut = new TProfile("fAvg_Cent_vs_Vy_Cent_woCut","<cent> vs Vy",100,VyCut[0], VyCut[1]);
799  fAvg_Cent_vs_Vy_Peri_woCut = new TProfile("fAvg_Cent_vs_Vy_Peri_woCut","<peri> vs Vy",100,VyCut[0], VyCut[1]);
801 
802  fAvg_Cent_vs_Vz_Cent_wCuts = new TProfile("fAvg_Cent_vs_Vz_Cent_wCuts","<cent> vs Vz",200,-10,10);
804  fAvg_Cent_vs_Vz_Peri_wCuts = new TProfile("fAvg_Cent_vs_Vz_Peri_wCuts","<peri> vs Vz",200,-10,10);
806  fAvg_Cent_vs_Vx_Cent_wCuts = new TProfile("fAvg_Cent_vs_Vx_Cent_wCuts","<cent> vs Vx",100,VxCut[0], VxCut[1]);
808  fAvg_Cent_vs_Vx_Peri_wCuts = new TProfile("fAvg_Cent_vs_Vx_Peri_wCuts","<peri> vs Vx",100,VxCut[0], VxCut[1]);
810  fAvg_Cent_vs_Vy_Cent_wCuts = new TProfile("fAvg_Cent_vs_Vy_Cent_wCuts","<cent> vs Vy",100,VyCut[0], VyCut[1]);
812  fAvg_Cent_vs_Vy_Peri_wCuts = new TProfile("fAvg_Cent_vs_Vy_Peri_wCuts","<peri> vs Vy",100,VyCut[0], VyCut[1]);
814 
815 
816  fTPCV0M_CentDiff_vs_Vz = new TProfile("fTPCV0M_CentDiff_vs_Vz","<TPC-V0M> vs Vz",200,-10,10);
818  fTPCV0M_CentDiff_vs_Vx = new TProfile("fTPCV0M_CentDiff_vs_Vx","<TPC-V0M> vs Vx",100,VxCut[0], VxCut[1]);
820  fTPCV0M_CentDiff_vs_Vy = new TProfile("fTPCV0M_CentDiff_vs_Vy","<TPC-V0M> vs Vy",100,VyCut[0], VyCut[1]);
822 
823  //std::cout<<"\n\n Rihan: AliAnalysisTaskVnZDC::CreateOutputObjects() is called....\n\n"<<std::endl;
824 
825  for(int i=0;i<frunflag;i++){
826  //store: ZDC energy for gain calibration:
827  fHist_ZDCA_En_Run[i] = new TProfile2D(Form("fHist_ZDCA_En_Run%d",runNums[i]),"",100,0,100,5,0,5,"");
828  fHist_ZDCC_En_Run[i] = new TProfile2D(Form("fHist_ZDCC_En_Run%d",runNums[i]),"",100,0,100,5,0,5,"");
829  //store: <X>,<Y> run by run for recenter://90 centrality bins:
830  for(int j=0;j<10;j++){
831  fHist_znCx_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znCx_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
832  fHist_znCy_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znCy_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
833  fHist_znAx_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znAx_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
834  fHist_znAy_V0_VxVy[i][j] = new TProfile2D(Form("fHist_znAy_V0_Run%d_Vz%d",runNums[i],j+1),"",NbinVt,0,NbinVt,90,0,90,"");
835  }
836  }
837 
838  for(int i=0;i<frunflag;i++){
841  for(int j=0;j<10;j++) {
842  fZDCESEList->Add(fHist_znCx_V0_VxVy[i][j]);
843  fZDCESEList->Add(fHist_znCy_V0_VxVy[i][j]);
844  fZDCESEList->Add(fHist_znAx_V0_VxVy[i][j]);
845  fZDCESEList->Add(fHist_znAy_V0_VxVy[i][j]);
846  }
847  }
848 
851 
852  }//---------- recenter1 ------------
853 
854 
855  fAnalysisUtil = new AliAnalysisUtils();
856  fAnalysisUtil->SetUseMVPlpSelection(kTRUE);
857  fAnalysisUtil->SetUseOutOfBunchPileUp(kTRUE);
858 
859 
860  PostData(1,fListHistos); //posting in slot #1,
861 
862  if(fAnalysisSet=="recenter1") {
863  PostData(2,fZDCESEList);
864  }
865  else{
866  fZListDummy = new TList();
867  fZListDummy->SetOwner(kTRUE);
870  PostData(2,fZListDummy);
871  }
872 
873  //AliDebug(2,
874  printf("\n\n::UserCreateOutPutObject(). NbinVt= %d, frunflag= %d, dataset: %s, Analysis= %s\n\n",NbinVt,frunflag,fDataSet.Data(),fAnalysisSet.Data());
875 
876 }
877 
878 //________________________________________________________________________
880 {
881 
882  int stepCount = 0;
883  //printf("\n ... ::UserExec() is being called. 1 Step %d... \n",stepCount);
884 
885  fHist_Event_count->Fill(stepCount);
886  stepCount++;
887 
888  /*if(!checkOnce){
889  checkOnce++;
890  }*/
891 
892  AliAODEvent *aod = dynamic_cast<AliAODEvent*>(InputEvent());
893  fEvent = dynamic_cast<AliFlowEventSimple*>(GetInputData(1));
894 
895 
896  if(!aod){
897  printf("\n ... ::UserExec no aod found..... \n");
898  return;
899  }
900 
901  fHist_Event_count->Fill(stepCount);
902  stepCount++;
903 
904 
905  AliAODVertex *pVertex = aod->GetPrimaryVertex();
906  Double_t Vxyz[3] = {0,0,0};
907  Vxyz[0] = pVertex->GetX();
908  Vxyz[1] = pVertex->GetY();
909  Vxyz[2] = pVertex->GetZ();
910 
911 
912  //------- Apply Necessary event cuts ---------
913  if(Vxyz[2] >= VzCut[1] || Vxyz[2] <= VzCut[0]) return;
914  fHist_Event_count->Fill(stepCount);
915  stepCount++;
916 
917  if(Vxyz[0] >= VxCut[1] || Vxyz[0] <= VxCut[0]) return;
918  fHist_Event_count->Fill(stepCount);
919  stepCount++;
920 
921  if(Vxyz[1] >= VyCut[1] || Vxyz[1] <= VyCut[0]) return;
922  fHist_Event_count->Fill(stepCount);
923  stepCount++;
924 
925 
926 
927  Double_t EvtCent = fEvent->GetCentrality();
928  fHist_Cent_count1-> Fill(EvtCent);
929 
930  //---------- get runindex: --------------
931  Int_t runNumber = aod->GetRunNumber();
932  Int_t runindex = -111;
933 
934  for(int i=0;i<frunflag;i++){
935  if(runNumber==runNums[i])
936  {
937  runindex = i;
938  break;
939  }
940  }
941  if(runindex<0) {
942  AliDebug(2,"\n ::UserExec Runnumber is not listed ..... \n");
943  exit(1);
944  }
945  //-----------------------------------------
946 
947 
948  //--------- starting pileup rejection work: --------
949  Double_t centrV0M=300; Double_t centrCL1=300;
950  Double_t centrCL0=300; Double_t centrTRK=300;
951 
952  if(fDataSet=="2010"||fDataSet=="2011"){
953  centrV0M = ((AliVAODHeader*)aod->GetHeader())->GetCentralityP()->GetCentralityPercentile("V0M");
954  centrCL1 = ((AliVAODHeader*)aod->GetHeader())->GetCentralityP()->GetCentralityPercentile("CL1");
955  centrCL0 = ((AliVAODHeader*)aod->GetHeader())->GetCentralityP()->GetCentralityPercentile("CL0");
956  centrTRK = ((AliVAODHeader*)aod->GetHeader())->GetCentralityP()->GetCentralityPercentile("TRK");
957  }// 2010/2011
958 
959  else{
960  fMultSelection = (AliMultSelection*) InputEvent()->FindListObject("MultSelection");
961  if(!fMultSelection) {
962  //If you get this warning, please check if AliMultSelectionTask actually ran (before your task)
963  AliDebug(2,Form("\n **WARNING** ::UserExec() AliMultSelection object not found. Step# %d\n",stepCount));
964  exit(1);
965  }
966  centrV0M = fMultSelection->GetMultiplicityPercentile("V0M");
967  centrCL1 = fMultSelection->GetMultiplicityPercentile("CL1");
968  centrCL0 = fMultSelection->GetMultiplicityPercentile("CL0");
969  centrTRK = fMultSelection->GetMultiplicityPercentile("TRK");
970  }// 2015
971 
972 
973  Bool_t BisPileup=kFALSE;
974 
975  if(fRejectPileUp && InputEvent()) {
976  //if(!fCutsEvent->IsSelected(InputEvent(),MCEvent())) return;
977  //if(fDataSet!=k2015 && fDataSet!=k2015v6) { //jacopo
978  if(fDataSet!="2015") {
979  if(plpMV(aod)) {
980  fPileUpCount->Fill(0.5);
981  BisPileup=kTRUE;
982  }
983  Int_t isPileup = aod->IsPileupFromSPD(3);
984  if(isPileup != 0) {
985  fPileUpCount->Fill(1.5);
986  BisPileup=kTRUE;
987  }
988  if(((AliAODHeader*)aod->GetHeader())->GetRefMultiplicityComb08() < 0) {
989  fPileUpCount->Fill(2.5);
990  BisPileup=kTRUE;
991  }
992  if(aod->IsIncompleteDAQ()) {
993  fPileUpCount->Fill(3.5);
994  BisPileup=kTRUE;
995  }
996 
997  //check vertex consistency
998  const AliAODVertex* vtTrc = aod->GetPrimaryVertex();
999  const AliAODVertex* vtSPD = aod->GetPrimaryVertexSPD();
1000 
1001  if(vtTrc->GetNContributors() < 2 || vtSPD->GetNContributors()<1) {
1002  fPileUpCount->Fill(5.5);
1003  BisPileup=kTRUE;
1004  }
1005 
1006  double covTrc[6], covSPD[6];
1007  vtTrc->GetCovarianceMatrix(covTrc);
1008  vtSPD->GetCovarianceMatrix(covSPD);
1009 
1010  double dz = vtTrc->GetZ() - vtSPD->GetZ();
1011 
1012  double errTot = TMath::Sqrt(covTrc[5]+covSPD[5]);
1013  double errTrc = TMath::Sqrt(covTrc[5]);
1014  double nsigTot = dz/errTot;
1015  double nsigTrc = dz/errTrc;
1016 
1017  if(TMath::Abs(dz)>0.2 || TMath::Abs(nsigTot)>10 || TMath::Abs(nsigTrc)>20) {
1018  fPileUpCount->Fill(6.5);
1019  BisPileup=kTRUE;
1020  }
1021  if(fAnalysisUtil->IsPileUpEvent(InputEvent())) {
1022  fPileUpCount->Fill(7.5);
1023  BisPileup=kTRUE;
1024  }
1025  }
1026 
1027  else {
1028  // pileup from AliMultSelection
1029  if(!fMultSelection->GetThisEventIsNotPileup())
1030  fPileUpMultSelCount->Fill(0.5);
1031  if(!fMultSelection->GetThisEventIsNotPileupMV())
1032  fPileUpMultSelCount->Fill(1.5);
1033  if(!fMultSelection->GetThisEventIsNotPileupInMultBins())
1034  fPileUpMultSelCount->Fill(2.5);
1035  if(!fMultSelection->GetThisEventHasNoInconsistentVertices())
1036  fPileUpMultSelCount->Fill(3.5);
1037  if(!fMultSelection->GetThisEventPassesTrackletVsCluster())
1038  fPileUpMultSelCount->Fill(4.5);
1039  if(!fMultSelection->GetThisEventIsNotAsymmetricInVZERO())
1040  fPileUpMultSelCount->Fill(5.5);
1041  if(!fMultSelection->GetThisEventIsNotIncompleteDAQ())
1042  fPileUpMultSelCount->Fill(6.5);
1043  if(!fMultSelection->GetThisEventHasGoodVertex2016())
1044  fPileUpMultSelCount->Fill(7.5);
1045 
1046  BisPileup=kFALSE;
1047 
1048  // pile-up a la Dobrin for LHC15o
1049  if(plpMV(aod)) {
1050  fPileUpCount->Fill(0.5);
1051  BisPileup=kTRUE;
1052  }
1053 
1054  Int_t isPileup = aod->IsPileupFromSPD(3);
1055  if(isPileup != 0) {
1056  fPileUpCount->Fill(1.5);
1057  BisPileup=kTRUE;
1058  }
1059 
1060  if(((AliAODHeader*)aod->GetHeader())->GetRefMultiplicityComb08() < 0) {
1061  fPileUpCount->Fill(2.5);
1062  BisPileup=kTRUE;
1063  }
1064 
1065  if(aod->IsIncompleteDAQ()) {
1066  fPileUpCount->Fill(3.5);
1067  BisPileup=kTRUE;
1068  }
1069 
1070  if(fabs(centrV0M-centrCL1)>7.5) {
1071  fPileUpCount->Fill(4.5);
1072  BisPileup=kTRUE;
1073  }
1074 
1075  // check vertex consistency
1076  const AliAODVertex* vtTrc = aod->GetPrimaryVertex();
1077  const AliAODVertex* vtSPD = aod->GetPrimaryVertexSPD();
1078 
1079  if (vtTrc->GetNContributors() < 2 || vtSPD->GetNContributors()<1) {
1080  fPileUpCount->Fill(5.5);
1081  BisPileup=kTRUE;
1082  }
1083 
1084  double covTrc[6], covSPD[6];
1085  vtTrc->GetCovarianceMatrix(covTrc);
1086  vtSPD->GetCovarianceMatrix(covSPD);
1087 
1088  double dz = vtTrc->GetZ() - vtSPD->GetZ();
1089 
1090  double errTot = TMath::Sqrt(covTrc[5]+covSPD[5]);
1091  double errTrc = TMath::Sqrt(covTrc[5]);
1092  double nsigTot = dz/errTot;
1093  double nsigTrc = dz/errTrc;
1094 
1095  if(TMath::Abs(dz)>0.2 || TMath::Abs(nsigTot)>10 || TMath::Abs(nsigTrc)>20) {
1096  fPileUpCount->Fill(6.5);
1097  BisPileup=kTRUE;
1098  }
1099 
1100  // cuts on tracks
1101  const Int_t nTracks = aod->GetNumberOfTracks();
1102  Int_t multEsd = ((AliAODHeader*)aod->GetHeader())->GetNumberOfESDTracks();
1103 
1104  Int_t multTrk = 0;
1105  Int_t multTrkBefC = 0;
1106  Int_t multTrkTOFBefC = 0;
1107  Int_t multTPC = 0;
1108 
1109  for(Int_t it = 0; it < nTracks; it++) {
1110  AliAODTrack* aodTrk = (AliAODTrack*)aod->GetTrack(it);
1111  if(!aodTrk) {
1112  delete aodTrk;
1113  continue;
1114  }
1115 // if(aodTrk->TestFilterBit(32)){
1116 // multTrkBefC++;
1117 // if(TMath::Abs(aodTrk->GetTOFsignalDz()) <= 10. && aodTrk->GetTOFsignal() >= 12000. && aodTrk->GetTOFsignal() <= 25000.)
1118 // multTrkTOFBefC++;
1119 // if((TMath::Abs(aodTrk->Eta()) < 0.8) && (aodTrk->GetTPCNcls() >= 70) && (aodTrk->Pt() >= 0.2) && (aodTrk->Pt() < 20.))
1120 // multTrk++;
1121 // }
1122  if(aodTrk->TestFilterBit(128))
1123  multTPC++;
1124  } // end of for (Int_t it = 0; it < nTracks; it++)
1125 
1126  Double_t multTPCn = multTPC;
1127  Double_t multEsdn = multEsd;
1128  Double_t multESDTPCDif = multEsdn - multTPCn*3.38;
1129 
1130  if(multESDTPCDif > (fRejectPileUpTight?700.:15000.)) {
1131  fPileUpCount->Fill(7.5);
1132  BisPileup=kTRUE;
1133  }
1134 
1135  if(fRejectPileUpTight) {
1136  if(BisPileup==kFALSE) {
1137  if(!fMultSelection->GetThisEventIsNotPileup()) BisPileup=kTRUE;
1138  if(!fMultSelection->GetThisEventIsNotPileupMV()) BisPileup=kTRUE;
1139  if(!fMultSelection->GetThisEventIsNotPileupInMultBins()) BisPileup=kTRUE;
1140  if(!fMultSelection->GetThisEventHasNoInconsistentVertices()) BisPileup=kTRUE;
1141  if(!fMultSelection->GetThisEventPassesTrackletVsCluster()) BisPileup=kTRUE;
1142  if(!fMultSelection->GetThisEventIsNotIncompleteDAQ()) BisPileup=kTRUE;
1143  if(!fMultSelection->GetThisEventHasGoodVertex2016()) BisPileup=kTRUE;
1144  if(BisPileup) fPileUpCount->Fill(8.5);
1145  }
1146  }
1147  }//-> 2015
1148  }
1149 
1150  //----------- pile up rejection done ---------
1151 
1152  //printf("\n ... ::UserExec() is being called. 4 Step %d... \n",stepCount);
1153 
1154  if(BisPileup){
1155  AliDebug(2,Form("\n ::UserExec Pileup event found, skipping.\n Dataset: %s \n",fDataSet.Data()));
1156  return;
1157  }
1158 
1159  fHist_Event_count->Fill(stepCount);
1160  stepCount++;
1161 
1162  if(fAnalysisSet=="recenter1") {
1163  fTPCV0M_CentDiff_vs_Vx->Fill(Vxyz[0],centrTRK-centrV0M);
1164  fTPCV0M_CentDiff_vs_Vy->Fill(Vxyz[1],centrTRK-centrV0M);
1165  fTPCV0M_CentDiff_vs_Vz->Fill(Vxyz[2],centrTRK-centrV0M);
1166 
1167  if(EvtCent<20){
1168  fAvg_Cent_vs_Vx_Cent_woCut->Fill(Vxyz[0],EvtCent);
1169  fAvg_Cent_vs_Vy_Cent_woCut->Fill(Vxyz[1],EvtCent);
1170  fAvg_Cent_vs_Vz_Cent_woCut->Fill(Vxyz[2],EvtCent);
1171  }
1172  if(EvtCent>50){
1173  fAvg_Cent_vs_Vx_Peri_woCut->Fill(Vxyz[0],EvtCent);
1174  fAvg_Cent_vs_Vy_Peri_woCut->Fill(Vxyz[1],EvtCent);
1175  fAvg_Cent_vs_Vz_Peri_woCut->Fill(Vxyz[2],EvtCent);
1176  }
1177  }
1178 
1179  // ********** fZDCgain alpha = 0.50 instead of 0.35 *********
1180  AliAODZDC *aodZDC = aod->GetZDCData();
1181  Float_t fZDCGainAlpha = 0.500;
1182  Float_t energyZNC = (Float_t) (aodZDC->GetZNCEnergy());
1183  Float_t energyZPC = (Float_t) (aodZDC->GetZPCEnergy());
1184  Float_t energyZNA = (Float_t) (aodZDC->GetZNAEnergy());
1185  Float_t energyZPA = (Float_t) (aodZDC->GetZPAEnergy());
1186  Float_t energyZEM1 = (Float_t) (aodZDC->GetZEM1Energy());
1187  Float_t energyZEM2 = (Float_t) (aodZDC->GetZEM2Energy());
1188 
1189  const Double_t * towZNC = aodZDC->GetZNCTowerEnergy();
1190  const Double_t * towZPC = aodZDC->GetZPCTowerEnergy();
1191  const Double_t * towZNA = aodZDC->GetZNATowerEnergy();
1192  const Double_t * towZPA = aodZDC->GetZPATowerEnergy();
1193 
1194  const Double_t * towZNClg = aodZDC->GetZNCTowerEnergyLR(); // Low gain something, should not be used.
1195  const Double_t * towZNAlg = aodZDC->GetZNATowerEnergyLR();
1196 
1197  Double_t towZPClg[5] = {0.,};
1198  Double_t towZPAlg[5] = {0.,};
1199 
1200  for(Int_t it=0; it<5; it++) {
1201  towZPClg[it] = 8*towZPC[it];
1202  towZPAlg[it] = 8*towZNA[it];
1203  }
1204 
1205  //sanity: remove if any of ZDC_C_A has negetive Energy:
1206  if(towZNC[1]<0 || towZNC[2]<0 || towZNC[3]<0 || towZNC[4]<0) return;
1207  fHist_Event_count->Fill(stepCount);
1208  stepCount++;
1209 
1210  if(towZNA[1]<0 || towZNA[2]<0 || towZNA[3]<0 || towZNA[4]<0) return;
1211  fHist_Event_count->Fill(stepCount);
1212  stepCount++;
1213 
1214  if(fAnalysisSet=="recenter1"){
1215  for(Int_t it=0; it<5; it++) {
1216  fHist_ZDCA_En_Run[runindex]->Fill(EvtCent,it,towZNA[it]);
1217  fHist_ZDCC_En_Run[runindex]->Fill(EvtCent,it,towZNC[it]);
1218  }
1219  }
1220 
1221 //********** Get centroid from ZDCs **************
1222 
1223  Double_t xyZNC[2]={999.,999.};
1224  Double_t xyZNA[2]={999.,999.};
1225 
1226  Float_t zncEnergy=0., znaEnergy=0.;
1227 
1228  for(Int_t i=0; i<5; i++){
1229  zncEnergy += towZNC[i];
1230  znaEnergy += towZNA[i];
1231  }
1232 
1233 
1234  Double_t AvTowerGain[8] = {1., 1., 1., 1., 1., 1., 1., 1.};
1235 
1236 /*---------------------------------------------------------------
1237  Int_t CenBin = GetCenBin(centrperc);
1238  Double_t zvtxpos[3]={0.,0.,0.};
1239  fFlowEvent->GetVertexPosition(zvtxpos);
1240  Int_t RunNum=fFlowEvent->GetRun();
1241  if(fTowerEqList) {
1242  if(RunNum!=fCachedRunNum) {
1243  for(Int_t i=0; i<8; i++) {
1244  fTowerGainEq[i] = (TH1D*)(fTowerEqList->FindObject(Form("Run %d",RunNum))->FindObject(Form("fhnTowerGainEqFactor[%d][%d]",RunNum,i)));
1245  }
1246  }
1247  }
1248  Bool_t fUseMCCen = kFALSE; //rihan:hardcoded
1249  if (fUseMCCen) {
1250  if(aod->GetRunNumber() < 209122)
1251  aodZDC->GetZNCentroidInPbPb(1380., xyZNC, xyZNA);
1252  else
1253  aodZDC->GetZNCentroidInPbPb(2510., xyZNC, xyZNA);
1254  }
1255  else {
1256  //set tower gain equalization, if available
1257  if(fTowerEqList) {
1258  for(Int_t i=0; i<8; i++)
1259  {
1260  if(fTowerGainEq[i])
1261  AvTowerGain[i] = fTowerGainEq[i]->GetBinContent(fTowerGainEq[i]->FindBin(centrperc));
1262  }
1263  }//--------------------------------------------------------------- */
1264 
1265 
1266 
1267  const Float_t x[4] = {-1.75, 1.75,-1.75, 1.75};
1268  const Float_t y[4] = {-1.75, -1.75, 1.75, 1.75};
1269 
1270  Float_t numXZNC=0., numYZNC=0., denZNC=0., wZNC;
1271  Float_t numXZNA=0., numYZNA=0., denZNA=0., wZNA;
1272 
1273  for(Int_t i=0; i<4; i++)
1274  {
1275  if(towZNC[i+1]>0.)
1276  {
1277  wZNC = TMath::Power(towZNC[i+1], fZDCGainAlpha)*AvTowerGain[i];
1278  numXZNC += x[i]*wZNC;
1279  numYZNC += y[i]*wZNC;
1280  denZNC += wZNC;
1281  }
1282 
1283  if(towZNA[i+1]>0.) {
1284  wZNA = TMath::Power(towZNA[i+1], fZDCGainAlpha)*AvTowerGain[i+4];
1285  numXZNA += x[i]*wZNA;
1286  numYZNA += y[i]*wZNA;
1287  denZNA += wZNA;
1288  }
1289  }
1290 
1291  if(denZNC!=0) {
1292  xyZNC[0] = numXZNC/denZNC;
1293  xyZNC[1] = numYZNC/denZNC;
1294  }
1295  else{
1296  xyZNC[0] = 999.;
1297  xyZNC[1] = 999.;
1298  zncEnergy = 0.;
1299  }
1300  if(denZNA!=0) {
1301  xyZNA[0] = numXZNA/denZNA;
1302  xyZNA[1] = numYZNA/denZNA;
1303  }
1304  else{
1305  xyZNA[0] = 999.;
1306  xyZNA[1] = 999.;
1307  znaEnergy = 0.;
1308  }
1309 
1310 
1311  //----- Important: zdcA_X = -zdcA_X ---------
1312  xyZNA[0] = -1.*xyZNA[0];
1313  //===========================================
1314 
1315 
1316  if(sqrt(xyZNC[0]*xyZNC[0] + xyZNC[1]*xyZNC[1])>1.5) return;
1317  fHist_Event_count->Fill(stepCount);
1318  stepCount++;
1319 
1320  if(sqrt(xyZNA[0]*xyZNA[0] + xyZNA[1]*xyZNA[1])>1.5) return;
1321  fHist_Event_count->Fill(stepCount);
1322  stepCount++;
1323 
1324  // ---------------- ZDC multiplicities posted at the ende of this code-----------------------
1325 
1326  fHist_EventperRun->Fill(runindex);
1327 
1328 
1329  //----- calculate RefMult for event ---------
1330  AliFlowTrackSimple* pTrack = NULL;
1331  Int_t iTracks = fEvent->NumberOfTracks();
1332 
1333  Double_t Qnx_TPC[3] = {0,};
1334  Double_t Qny_TPC[3] = {0,};
1335  Double_t psi2,dPhi,dPt,dEta;
1336  Int_t ipTBin;
1337  Double_t pTwgt;
1338  Int_t nRefMult = 0;
1339  Double_t npoiMult = 0;
1340 
1341  Int_t cIndex = -1;
1342  if(EvtCent<5.0){ cIndex = 0;}
1343  else if(EvtCent>=5.0 && EvtCent<10){
1344  cIndex = 1;
1345  }
1346  else {
1347  cIndex = abs(EvtCent/10.0) + 1;
1348  }
1349 
1350 
1351  Double_t dUx = 0;
1352  Double_t dUy = 0;
1353 
1354  if(fAnalysisSet=="recenter2") {
1355  for(Int_t i=0; i<iTracks; i++)
1356  {
1357  pTrack = fEvent->GetTrack(i);
1358  if (!pTrack) continue;
1359  dPhi = pTrack->Phi();
1360  dPt = pTrack-> Pt();
1361  dEta = pTrack->Eta();
1362  //dCharge = pTrack->Charge();
1363  if(fabs(dEta)>0.8) continue;
1364  if(dPt<0.15 || dPt>10.0) continue;
1365  nRefMult++;
1366  if(dPt<0.20) continue;
1367 
1368  //fHistQA_etaphi->Fill(dPhi,dEta);
1369  ipTBin = fFB_Efficiency_Cent[cIndex]->FindBin(dPt);
1370  pTwgt = 1.0/fFB_Efficiency_Cent[cIndex]->GetBinContent(ipTBin);
1371 
1372  Qnx_TPC[0] += TMath::Cos(2.*dPhi)*pTwgt;
1373  Qny_TPC[0] += TMath::Sin(2.*dPhi)*pTwgt;
1374 
1375  npoiMult += pTwgt;
1376  //cout<<"c = "<<EvtCent<<" indC = "<<cIndex<<" "<<npoiMult<<" pt = "<<dPt<<" weight = "<<pTwgt<<endl;
1377  }//track loop ends
1378  } //----- call track loop only for recenter2 -------
1379 
1380 
1381  if(npoiMult>0){
1382  dUx = Qnx_TPC[0]/npoiMult;
1383  dUy = Qny_TPC[0]/npoiMult;
1384  }
1385  else{
1386  dUx = 0;
1387  dUy = 0;
1388  }
1389 
1390  //double psi2 = 0.5*TMath::ATan2(Qny_TPC[0],Qnx_TPC[0]);
1391  //if(psi2<0) psi2 += TMath::Pi();
1392  //fHist_EventPlane2->Fill(psi2);
1393  //fill q vectors for TPC recentering
1394  //fHist_QnxRecent->Fill(EvtCent,runindex,(Qnx_TPC[0]/(double)npoiMult));
1395  //fHist_QnyRecent->Fill(EvtCent,runindex,(Qny_TPC[0]/(double)npoiMult));
1396 
1397  fHist_Cent_count2-> Fill(EvtCent);
1398 
1399 
1400  if(fAnalysisSet=="recenter1") {
1401  if(EvtCent<20){
1402  fAvg_Cent_vs_Vx_Cent_wCuts->Fill(Vxyz[0],EvtCent);
1403  fAvg_Cent_vs_Vy_Cent_wCuts->Fill(Vxyz[1],EvtCent);
1404  fAvg_Cent_vs_Vz_Cent_wCuts->Fill(Vxyz[2],EvtCent);
1405  }
1406  if(EvtCent>50){
1407  fAvg_Cent_vs_Vx_Peri_wCuts->Fill(Vxyz[0],EvtCent);
1408  fAvg_Cent_vs_Vy_Peri_wCuts->Fill(Vxyz[1],EvtCent);
1409  fAvg_Cent_vs_Vz_Peri_wCuts->Fill(Vxyz[2],EvtCent);
1410  }
1411  }
1412 
1413  Int_t nTracks = aod->GetNumberOfTracks(); //number of AOD tracks
1414 
1415  fHist_Vertex_Vz->Fill(Vxyz[2]);
1416  fHist_Vertex_XY->Fill(Vxyz[0],Vxyz[1]);
1417 
1418  fHist_Vx_vs_runnum ->Fill(runindex,Vxyz[0]);
1419  fHist_Vy_vs_runnum ->Fill(runindex,Vxyz[1]);
1420  fHist_Vz_vs_runnum ->Fill(runindex,Vxyz[2]);
1421  fHist_tracks_vs_runnum->Fill(runindex,nTracks);
1422 
1423 
1424 
1425 
1426  psi2 = TMath::ATan2(xyZNC[1],xyZNC[0]);
1427  if(psi2<0) psi2 += 2.*TMath::Pi();
1428  fHist_Psi1_zdnC->Fill(psi2);
1429 
1430  psi2 = TMath::ATan2(xyZNA[1],xyZNA[0]);
1431  if(psi2<0) psi2 += 2.*TMath::Pi();
1432  fHist_Psi1_zdnA->Fill(psi2);
1433 
1434 
1435 
1436  Int_t indexVx = fHist_Vx_ArrayFinder->FindBin(Vxyz[0]);
1437  Int_t indexVy = fHist_Vy_ArrayFinder->FindBin(Vxyz[1]);
1438  Int_t indexVz = fHist_Vz_ArrayFinder->FindBin(Vxyz[2]);
1439 
1440 
1441  /*if(fAnalysisSet=="recenter1"){
1442  Double_t tVertexBin1 = (Double_t) (indexVy-1)*vxBin + (Double_t)indexVx - 0.5 ; //
1443  fHist_vBincount->Fill(tVertexBin1);
1444  fHist_znCx_V0_VxVy[runindex][indexVz-1]->Fill(tVertexBin1,EvtCent,xyZNC[0]); //EvtCent
1445  fHist_znCy_V0_VxVy[runindex][indexVz-1]->Fill(tVertexBin1,EvtCent,xyZNC[1]);
1446  fHist_znAx_V0_VxVy[runindex][indexVz-1]->Fill(tVertexBin1,EvtCent,xyZNA[0]);
1447  fHist_znAy_V0_VxVy[runindex][indexVz-1]->Fill(tVertexBin1,EvtCent,xyZNA[1]);
1448  } */
1449 
1450  Int_t tVertexBin2 = (indexVy-1)*vxBin + indexVx; // 'Int_t' because bin starts with 1.
1451  Int_t tCentBin = abs(EvtCent) + 1;
1452  Double_t CentWgt = 1.0;
1453 
1454 
1455  //printf("\n ... ::UserExec() is being called. Step 6 %d... \n",stepCount);
1456 
1457  Double_t meanCx[2] = {0,};
1458  Double_t meanCy[2] = {0,};
1459  Double_t meanAx[2] = {0,};
1460  Double_t meanAy[2] = {0,};
1461 
1462  if(fAnalysisSet=="recenter2") {
1463 
1464  Double_t tVertexBinD = (Double_t)tVertexBin2 - 0.5;
1465  fHist_vBincount->Fill(tVertexBinD);
1466 
1467  Double_t FillVsWith[5] = {EvtCent,static_cast<Double_t>(nRefMult), Vxyz[0], Vxyz[1], Vxyz[2]};
1468  Double_t FillValue1[4] = {xyZNA[0],xyZNC[0],xyZNA[1],xyZNC[1]};
1469  Double_t FillValue11[4] = {xyZNA[0]*xyZNC[0],xyZNA[1]*xyZNC[1],xyZNC[0]*xyZNA[1],xyZNC[1]*xyZNA[0]}; //XaXc,YaYc,XcYa,YcXa
1470 
1471  //fill the uncorrected Qns first:
1472  for(int i=0;i<4;i++){
1473  for(int j=0;j<5;j++){
1474  fHist_X_vs_Obs_before[i][j]->Fill(FillVsWith[j],FillValue1[i]);
1475  fHist_XX_vs_Obs_before[i][j]->Fill(FillVsWith[j],FillValue11[i]);
1476  }
1477  }
1478 
1479  meanCx[0] = fHist_znCx_V0_VxVy[runindex][indexVz-1]->GetBinContent(tVertexBin2,tCentBin);
1480  meanCy[0] = fHist_znCy_V0_VxVy[runindex][indexVz-1]->GetBinContent(tVertexBin2,tCentBin);
1481  meanAx[0] = fHist_znAx_V0_VxVy[runindex][indexVz-1]->GetBinContent(tVertexBin2,tCentBin);
1482  meanAy[0] = fHist_znAy_V0_VxVy[runindex][indexVz-1]->GetBinContent(tVertexBin2,tCentBin);
1483 
1484 
1485 
1486 
1487 
1488  xyZNC[0] = xyZNC[0] - meanCx[0];
1489  xyZNC[1] = xyZNC[1] - meanCy[0];
1490  xyZNA[0] = xyZNA[0] - meanAx[0];
1491  xyZNA[1] = xyZNA[1] - meanAy[0];
1492 
1493  Double_t FillValue2[4] = {xyZNA[0],xyZNC[0],xyZNA[1],xyZNC[1]};
1494  Double_t FillValue21[4] = {xyZNA[0]*xyZNC[0],xyZNA[1]*xyZNC[1],xyZNC[0]*xyZNA[1],xyZNC[1]*xyZNA[0]}; //XaXc,YaYc,XcYa,YcXa
1495 
1496  //fill QA after recenter 1:
1497  for(int i=0;i<4;i++){
1498  for(int j=0;j<5;j++){
1499  fHist_X_vs_Obs_after1[i][j]->Fill(FillVsWith[j],FillValue2[i]);
1500  fHist_XX_vs_Obs_after1[i][j]->Fill(FillVsWith[j],FillValue21[i]);
1501  }
1502  }
1503 
1504  psi2 = TMath::ATan2(xyZNC[1],xyZNC[0]);
1505  if(psi2<0) psi2 += 2.*TMath::Pi();
1506  fHist_Psi1_zdnC_after1->Fill(psi2);
1507 
1508  psi2 = TMath::ATan2(xyZNA[1],xyZNA[0]);
1509  if(psi2<0) psi2 += 2.*TMath::Pi();
1510  fHist_Psi1_zdnA_after1->Fill(psi2);
1511 
1512  //fill results histograms:
1513  fHist_v2xV1_ZDN_Norm_cosXX ->Fill(EvtCent, dUx*xyZNA[0]*xyZNC[0],CentWgt);
1514  fHist_v2xV1_ZDN_Cent_cosXX ->Fill(EvtCent, dUx*xyZNA[0]*xyZNC[0],CentWgt);
1515  fHist_v2xV1_ZDN_Refm_cosXX ->Fill(nRefMult,dUx*xyZNA[0]*xyZNC[0],CentWgt);
1516 
1517  fHist_v2xV1_ZDN_Norm_cosYY ->Fill(EvtCent, dUx*xyZNA[1]*xyZNC[1],CentWgt);
1518  fHist_v2xV1_ZDN_Cent_cosYY ->Fill(EvtCent, dUx*xyZNA[1]*xyZNC[1],CentWgt);
1519  fHist_v2xV1_ZDN_Refm_cosYY ->Fill(nRefMult,dUx*xyZNA[1]*xyZNC[1],CentWgt);
1520 
1521  fHist_v2xV1_ZDN_Norm_sinXY ->Fill(EvtCent, dUy*xyZNA[0]*xyZNC[1],CentWgt);
1522  fHist_v2xV1_ZDN_Cent_sinXY ->Fill(EvtCent, dUy*xyZNA[0]*xyZNC[1],CentWgt);
1523  fHist_v2xV1_ZDN_Refm_sinXY ->Fill(nRefMult,dUy*xyZNA[0]*xyZNC[1],CentWgt);
1524 
1525  fHist_v2xV1_ZDN_Norm_sinYX ->Fill(EvtCent, dUy*xyZNA[1]*xyZNC[0],CentWgt);
1526  fHist_v2xV1_ZDN_Cent_sinYX ->Fill(EvtCent, dUy*xyZNA[1]*xyZNC[0],CentWgt);
1527  fHist_v2xV1_ZDN_Refm_sinYX ->Fill(nRefMult,dUy*xyZNA[1]*xyZNC[0],CentWgt);
1528 
1529  fHist_ZDN_resol_Norm_cos ->Fill(EvtCent, xyZNA[0]*xyZNC[0],CentWgt);
1530  fHist_ZDN_resol_Cent_cos ->Fill(EvtCent, xyZNA[0]*xyZNC[0],CentWgt);
1531  fHist_ZDN_resol_Refm_cos ->Fill(nRefMult,xyZNA[0]*xyZNC[0],CentWgt);
1532 
1533  fHist_ZDN_resol_Norm_sin ->Fill(EvtCent, xyZNA[1]*xyZNC[1],CentWgt);
1534  fHist_ZDN_resol_Cent_sin ->Fill(EvtCent, xyZNA[1]*xyZNC[1],CentWgt);
1535  fHist_ZDN_resol_Refm_sin ->Fill(nRefMult,xyZNA[1]*xyZNC[1],CentWgt);
1536 
1537  Double_t fullTerm = dUx*xyZNA[0]*xyZNC[0]-dUx*xyZNA[1]*xyZNC[1]+dUy*xyZNA[0]*xyZNC[1]+dUy*xyZNA[1]*xyZNC[0];
1538  Double_t fullReso = xyZNA[0]*xyZNC[0]+xyZNA[1]*xyZNC[1];
1539 
1540  fHist_v2xV1_ZDN_Norm_All ->Fill(EvtCent,fullTerm,CentWgt);
1541  fHist_v2xV1_ZDN_Cent_All ->Fill(EvtCent,fullTerm,CentWgt);
1542  fHist_v2xV1_ZDN_Refm_All ->Fill(nRefMult,fullTerm,CentWgt);
1543  fHist_ZDN_resol_Norm_All ->Fill(EvtCent,fullReso,CentWgt);
1544  fHist_ZDN_resol_Cent_All ->Fill(EvtCent,fullReso,CentWgt);
1545  fHist_ZDN_resol_Refm_All ->Fill(nRefMult,fullReso,CentWgt);
1546  }
1547 
1548 
1549 
1550  //fHist_ZDCn_A_XYvsRun->Fill(xyZNA[0],xyZNA[1],runindex);
1551  //fHist_ZDCn_C_XYvsRun->Fill(xyZNC[0],xyZNC[1],runindex);
1552 
1553 
1554 
1555  //if(fievent%10==0){
1556  //std::cout<<fievent<<" cTPC= "<<EvtCent<<" cV0M = "<<centrV0M<<" cWgt = "<<CentWgt<<" vz= "<<Vxyz[2]<<"\tCx= "<<meanCx[0]
1557  //<<"\tCy= "<<meanCy[0]<<"\tAx= "<<meanAx[0]<<"\tAy= "<<meanAy[0]<<"\tRefm= "<<nRefMult<<std::endl; }
1558 
1559 
1560  PostData(1,fListHistos);
1561 
1562  if(fAnalysisSet=="recenter2"){
1563  PostData(2,fZListDummy);
1564  }
1565  else{
1566  PostData(2,fZDCESEList);
1567  }
1568 
1569 
1570 
1571  fievent++;
1572 
1573 }
1574 
1575 
1576 
1577 
1579 {
1580  // Called once at the end of the query
1581  /* AliFlowAnalysisIDCSP *fSPTerm = new AliFlowAnalysisIDCSP();
1582  fListHistos = (TList*) GetOutputData(1);
1583  if(fListHistos)
1584  {
1585  fSPTerm->GetOutputHistograms(fListHistos);
1586  fSPTerm->Finish();
1587  PostData(1,fListHistos);
1588  }
1589  else
1590  {
1591  std::cout << "histgram list pointer is empty in Scalar Product" << endl;
1592  } */
1593  AliDebug(2,"\n ... AliAnalysisTaskVnZDC::Terminate() is being called ... \n");
1594 }
1595 
1596 
1597 double AliAnalysisTaskVnZDC::GetWDist(const AliVVertex* v0, const AliVVertex* v1)
1598 {
1599  // calculate sqrt of weighted distance to other vertex
1600  if(!v0 || !v1) {
1601  AliDebug(2,"\n\n ::GetWDist => One of vertices is not valid\n\n");
1602  return 0;
1603  }
1604  static TMatrixDSym vVb(3);
1605  double dist = -1;
1606  double dx = v0->GetX()-v1->GetX();
1607  double dy = v0->GetY()-v1->GetY();
1608  double dz = v0->GetZ()-v1->GetZ();
1609  double cov0[6],cov1[6];
1610  v0->GetCovarianceMatrix(cov0);
1611  v1->GetCovarianceMatrix(cov1);
1612  vVb(0,0) = cov0[0]+cov1[0];
1613  vVb(1,1) = cov0[2]+cov1[2];
1614  vVb(2,2) = cov0[5]+cov1[5];
1615  vVb(1,0) = vVb(0,1) = cov0[1]+cov1[1];
1616  vVb(0,2) = vVb(1,2) = vVb(2,0) = vVb(2,1) = 0.;
1617  vVb.InvertFast();
1618  if (!vVb.IsValid()) {
1619  AliDebug(2,"Singular Matrix\n");
1620  return dist;}
1621  dist = vVb(0,0)*dx*dx + vVb(1,1)*dy*dy + vVb(2,2)*dz*dz
1622  + 2*vVb(0,1)*dx*dy + 2*vVb(0,2)*dx*dz + 2*vVb(1,2)*dy*dz;
1623  return dist>0 ? TMath::Sqrt(dist) : -1;
1624 }
1625 
1627  { // check for multi-vertexer pile-up
1628  const int kMinPlpContrib = 5;
1629  const double kMaxPlpChi2 = 5.0;
1630  const double kMinWDist = 15;
1631 
1632  const AliVVertex* vtPrm = 0;
1633  const AliVVertex* vtPlp = 0;
1634 
1635  int nPlp = 0;
1636 
1637  if(!(nPlp=aod->GetNumberOfPileupVerticesTracks()))
1638  return kFALSE;
1639 
1640  vtPrm = aod->GetPrimaryVertex();
1641  if(vtPrm == aod->GetPrimaryVertexSPD())
1642  return kTRUE; // there are pile-up vertices but no primary
1643 
1644  //int bcPrim = vtPrm->GetBC();
1645 
1646  for(int ipl=0;ipl<nPlp;ipl++) {
1647  vtPlp = (const AliVVertex*)aod->GetPileupVertexTracks(ipl);
1648  if(vtPlp->GetNContributors() < kMinPlpContrib) continue;
1649  if(vtPlp->GetChi2perNDF() > kMaxPlpChi2) continue;
1650  //int bcPlp = vtPlp->GetBC();
1651  //if (bcPlp!=AliVTrack::kTOFBCNA && TMath::Abs(bcPlp-bcPrim)>2)
1652  // return kTRUE; // pile-up from other BC
1653 
1654  double wDst = GetWDist(vtPrm,vtPlp);
1655  if(wDst<kMinWDist) continue;
1656 
1657  return kTRUE; // pile-up: well separated vertices
1658  }
1659 
1660  return kFALSE;
1661  }
1662 
1663 
1664 
1665 
ClassImp(AliAnalysisTaskVnZDC) AliAnalysisTaskVnZDC
TProfile * fHist_X_vs_Obs_before[4][5]
AliFlowEventSimple * fEvent
double Double_t
Definition: External.C:58
Definition: External.C:260
Definition: External.C:236
TProfile2D * fHist_znAx_V0_VxVy[90][10]
AliMultSelection * fMultSelection
AliFlowTrackSimple * GetTrack(Int_t i)
TProfile2D * fHist_znCy_V0_VxVy[90][10]
TProfile * fHist_XX_vs_Obs_before[4][5]
virtual void Terminate(Option_t *)
Int_t fievent
input event
Double_t Phi() const
int Int_t
Definition: External.C:63
float Float_t
Definition: External.C:68
Bool_t plpMV(const AliAODEvent *aod)
TProfile * fHist_XX_vs_Obs_after1[4][5]
TProfile2D * fHist_ZDCC_En_Run[90]
Definition: External.C:212
virtual void UserCreateOutputObjects()
Double_t GetCentrality() const
AliAnalysisUtils * fAnalysisUtil
TProfile * fHist_X_vs_Obs_after1[4][5]
double GetWDist(const AliVVertex *v0, const AliVVertex *v1)
TProfile2D * fHist_ZDCA_En_Run[90]
const char Option_t
Definition: External.C:48
TProfile2D * fHist_znAy_V0_VxVy[90][10]
TProfile2D * fHist_znCx_V0_VxVy[90][10]
Double_t Eta() const
bool Bool_t
Definition: External.C:53
virtual void UserExec(Option_t *option)
Int_t NumberOfTracks() const