AliPhysics  9b6b435 (9b6b435)
AliAnalysisTaskCMEV0PID.cxx
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  * *
4  * Author: The ALICE Off-line Project. *
5  * Contributors are mentioned in the code where appropriate. *
6  * *
7  * Permission to use, copy, modify and distribute this software and its *
8  * documentation strictly for non-commercial purposes is hereby granted *
9  * without fee, provided that the above copyright notice appears in all *
10  * copies and that both the copyright notice and this permission notice *
11  * appear in the supporting documentation. The authors make no claims *
12  * about the suitability of this software for any purpose. It is *
13  * provided "as is" without express or implied warranty. *
14  **************************************************************************/
15 
16 /* $Id: AliAnalysisTaskCMEV0PID.cxx Rihan Haque 14/02/2018 $ */
17 
18 //-- general include---
19 #include "TChain.h"
20 #include "TTree.h"
21 #include "TGrid.h"
22 #include "TROOT.h"
23 #include "TObjArray.h"
24 #include "TMatrixDSym.h"
25 
26 #include "TMath.h"
27 #include "stdio.h"
28 #include "Riostream.h"
29 
30 //---- manager and handler---
31 #include "AliAnalysisManager.h"
32 #include "AliInputEventHandler.h"
33 
34 //---V0 and ZDC info---
35 #include "AliAODZDC.h"
36 #include "AliAODVZERO.h"
37 #include "AliAODVertex.h"
38 
39 //---AOD,ESD event--
40 #include "AliESDEvent.h"
41 #include "AliAODHeader.h"
42 #include "AliAODEvent.h"
43 #include "AliAODTrack.h"
44 
45 //----for PID-----
46 #include "AliPIDResponse.h"
47 #include "AliPIDCombined.h"
48 
49 //----- Vevent and tracks
50 #include "AliVEventHandler.h"
51 #include "AliVEvent.h"
52 #include "AliVTrack.h"
53 #include "AliVParticle.h"
54 #include "AliCentrality.h"
55 
56 //----- must include-------
57 #include "AliMultSelection.h"
58 #include "AliAnalysisUtils.h"
59 #include "AliPhysicsSelection.h"
60 #include "AliFlowEventSimple.h"
61 #include "AliAnalysisTaskSE.h"
63 
64 //using namespace std;
65 
66 using std::cout;
67 using std::endl;
68 using std::vector;
69 
70 
72 
73 
75  fVevent(NULL),
76  fESD(NULL),
77  fAOD(NULL),
78  fPIDResponse(NULL),
79  fMultSelection(NULL),
80  fAnalysisUtil(NULL),
81  fListHist(NULL),
82  mfileFBHijing(NULL),
83  fListFBHijing(NULL),
84  fListNUACorr(NULL),
85  fListV0MCorr(NULL),
86  fHistTaskConfigParameters(NULL),
87  fHistPileUpCount(NULL),
88  fHistMultSelPUCount(NULL),
89  fHistEtaPtBefore(NULL),
90  fHistEtaPtAfter(NULL),
91  fHistTPCvsGlobalMultBefore(NULL),
92  fHistTPCvsGlobalMultAfter(NULL),
93  fHistTPCdEdxvsPBefore(NULL),
94  fHistTPCdEdxvsPAfter(NULL),
95  fHistTOFBetavsPBefore(NULL),
96  fHistTOFBetavsPAfter(NULL),
97  fHistTOFMassvsPtBefore(NULL),
98  fHistTOFMatchCount(NULL),
99  fHistTPCVsESDTrkBefore(NULL),
100  fHistTPCVsESDTrkAfter(NULL),
101  fHistTPConlyVsCL1Before(NULL),
102  fHistTPConlyVsV0MBefore(NULL),
103  fHistTPConlyVsCL1After(NULL),
104  fHistTPConlyVsV0MAfter(NULL),
105  fHistGlobalVsV0MBefore(NULL),
106  fHistGlobalVsV0MAfter(NULL),
107  fHistRawVsCorrMultFB(NULL),
108  hCentvsTPCmultCuts(NULL),
109  fHV0AEventPlaneVsCent(NULL),
110  fHV0CEventPlaneVsCent(NULL),
111  fHTPCAEventPlaneVsCent(NULL),
112  fHTPCCEventPlaneVsCent(NULL),
113  fHTPCEventPlaneVsCent(NULL),
114  fV0MultChVsRun(NULL),
115  fCentDistBefore(NULL),
116  fCentDistAfter(NULL),
117  fHCorrectV0M(NULL),
118  fHAvgerageQnV0A(NULL),
119  fHAvgerageQnV0C(NULL),
120  fHCentWeightForRun(NULL),
121  fQAEtaPhiAfterNUA(NULL),
122  fQAEtaPhiAfterNUAPion(NULL),
123  fQAEtaPhiAfterNUAKaon(NULL),
124  fQAEtaPhiAfterNUAProton(NULL),
125  fV0AQ2xVsCentRun(NULL),
126  fV0AQ2yVsCentRun(NULL),
127  fV0CQ2xVsCentRun(NULL),
128  fV0CQ2yVsCentRun(NULL),
129  fV0AQ3xVsCentRun(NULL),
130  fV0AQ3yVsCentRun(NULL),
131  fV0CQ3xVsCentRun(NULL),
132  fV0CQ3yVsCentRun(NULL),
133  fTPCAQ2xVsCentRun(NULL),
134  fTPCAQ2yVsCentRun(NULL),
135  fTPCCQ2xVsCentRun(NULL),
136  fTPCCQ2yVsCentRun(NULL),
137  fTPCAQ3xVsCentRun(NULL),
138  fTPCAQ3yVsCentRun(NULL),
139  fTPCCQ3xVsCentRun(NULL),
140  fTPCCQ3yVsCentRun(NULL),
141  fTPCAQ4xVsCentRun(NULL),
142  fTPCAQ4yVsCentRun(NULL),
143  fTPCCQ4xVsCentRun(NULL),
144  fTPCCQ4yVsCentRun(NULL),
145  fTPCFQ2xVsCentRun(NULL),
146  fTPCFQ2yVsCentRun(NULL),
147  fZPASignalPerChVsCent(NULL),
148  fZPCSignalPerChVsCent(NULL),
149  fZNASignalPerChVsCent(NULL),
150  fZNCSignalPerChVsCent(NULL),
151  fCentDistVsVz(NULL),
152  fFilterBit(1),
153  gN(1),
154  gM(1),
155  gPsiN(2),
156  fOldRunNum(111),
157  fEventCount(0),
158  fTPCclustMin(70),
159  fNSigmaCut(2),
160  fMinPtCut(0.2),
161  fMaxPtCut(5.0),
162  fMinEtaCut(-0.8),
163  fMaxEtaCut(0.8),
164  fTrkChi2Min(0.1),
165  fDCAxyMax(2.4),
166  fDCAzMax(3.2),
167  fdEdxMin(10.),
168  fCentralityPercentMin(0),
169  fCentralityPercentMax(90),
170  fPileUpSlopeParm(3.43),
171  fPileUpConstParm(43),
172  fMinVzCut(-10.0),
173  fMaxVzCut(10.0),
174  bApplyMCcorr(kFALSE),
175  bV0MGainCorr(kFALSE),
176  bSkipPileUpCut(kFALSE),
177  bFillNUAHistPID(kFALSE),
178  bUseKinkTracks(kFALSE),
179  sPathOfMCFile("/alien"),
180  sNucleiTP("PbPb"),
181  sCentrEstimator("V0M"),
182  fHistEventCount(NULL)
183 {
184  for(int i=0;i<3;i++){
185  fHistPtwithTPCNsigma[i]=NULL;
186  fHistPtwithTOFmasscut[i]=NULL;
187  fHistPtwithTOFSignal[i]=NULL;
188  fHistTOFnSigmavsPtAfter[i]=NULL;
189  fHistTPCnSigmavsPtAfter[i]=NULL;
190  fHistTPCTOFnSigmavsPtAfter[i]=NULL;
191  fHistTPCdEdxvsPtPIDAfter[i]=NULL;
192  }
193  for(int i=0;i<5;i++){
194  fHCorrectNUApos[i] = NULL;
195  fHCorrectNUAneg[i] = NULL;
196  }
197  for(int i=0;i<5;i++){ // for PID
198  fHCorrectNUAposPion[i] = NULL;
199  fHCorrectNUAnegPion[i] = NULL;
200  fHCorrectNUAposKaon[i] = NULL;
201  fHCorrectNUAnegKaon[i] = NULL;
202  fHCorrectNUAposProton[i] = NULL;
203  fHCorrectNUAnegProton[i] = NULL;
204  }
205  //3p vs centrality
206  for(int i=0;i<2;i++){
207  for(int j=0;j<4;j++){
208  fHist_Corr3p_EP_Norm_PN[i][j] = NULL;
209  fHist_Corr3p_EP_Norm_PP[i][j] = NULL;
210  fHist_Corr3p_EP_Norm_NN[i][j] = NULL;
211  }
212  for(int j=0;j<4;j++) {
213  fHist_Reso2n_EP_Norm_Det[i][j] = NULL;
214  }
215  }
216 
217  //Differential 3p Charge:
218  for(int i=0;i<2;i++){
219  for(int j=0;j<6;j++){
220  fHist_Corr3p_pTSum_EP_V0A_PN[i][j] = NULL;
221  fHist_Corr3p_pTSum_EP_V0A_PP[i][j] = NULL;
222  fHist_Corr3p_pTSum_EP_V0A_NN[i][j] = NULL;
223  fHist_Corr3p_pTSum_EP_V0C_PN[i][j] = NULL;
224  fHist_Corr3p_pTSum_EP_V0C_PP[i][j] = NULL;
225  fHist_Corr3p_pTSum_EP_V0C_NN[i][j] = NULL;
226 
227  fHist_Corr3p_pTDiff_EP_V0A_PN[i][j] = NULL;
228  fHist_Corr3p_pTDiff_EP_V0A_PP[i][j] = NULL;
229  fHist_Corr3p_pTDiff_EP_V0A_NN[i][j] = NULL;
230  fHist_Corr3p_pTDiff_EP_V0C_PN[i][j] = NULL;
231  fHist_Corr3p_pTDiff_EP_V0C_PP[i][j] = NULL;
232  fHist_Corr3p_pTDiff_EP_V0C_NN[i][j] = NULL;
233 
234  fHist_Corr3p_EtaDiff_EP_V0A_PN[i][j] = NULL;
235  fHist_Corr3p_EtaDiff_EP_V0A_PP[i][j] = NULL;
236  fHist_Corr3p_EtaDiff_EP_V0A_NN[i][j] = NULL;
237  fHist_Corr3p_EtaDiff_EP_V0C_PN[i][j] = NULL;
238  fHist_Corr3p_EtaDiff_EP_V0C_PP[i][j] = NULL;
239  fHist_Corr3p_EtaDiff_EP_V0C_NN[i][j] = NULL;
240  }
241  }
242 
243  for(int i=0;i<4;i++){
244  for(int j=0;j<5;j++){
245  fHist3DEtaPhiVz_Pos_Run[i][j]=NULL;
246  fHist3DEtaPhiVz_Neg_Run[i][j]=NULL;
247  }
248  }
249  for(int i=0;i<10;i++){
250  fFB_Efficiency_Cent[i] = NULL;
251  //fFB_Efficiency_Pion_Cent[i] = NULL;
252  //fFB_Efficiency_Kaon_Cent[i] = NULL;
253  //fFB_Efficiency_Proton_Pos_Cent[i] = NULL;
254  //fFB_Efficiency_Proton_Neg_Cent[i] = NULL;
255  }
256 
257 
258  DefineInput(0,TChain::Class());
259  DefineOutput(1,TList::Class());
260 }
261 //______________________________empty constructor_______________________
264  fVevent(NULL),
265  fESD(NULL),
266  fAOD(NULL),
267  fPIDResponse(NULL),
268  fMultSelection(NULL),
269  fAnalysisUtil(NULL),
270  fListHist(NULL),
271  mfileFBHijing(NULL),
272  fListFBHijing(NULL),
273  fListNUACorr(NULL),
274  fListV0MCorr(NULL),
275  fHistTaskConfigParameters(NULL),
276  fHistPileUpCount(NULL),
277  fHistMultSelPUCount(NULL),
278  fHistEtaPtBefore(NULL),
279  fHistEtaPtAfter(NULL),
280  fHistTPCvsGlobalMultBefore(NULL),
281  fHistTPCvsGlobalMultAfter(NULL),
282  fHistTPCdEdxvsPBefore(NULL),
283  fHistTPCdEdxvsPAfter(NULL),
284  fHistTOFBetavsPBefore(NULL),
285  fHistTOFBetavsPAfter(NULL),
286  fHistTOFMassvsPtBefore(NULL),
287  fHistTOFMatchCount(NULL),
288  fHistTPCVsESDTrkBefore(NULL),
289  fHistTPCVsESDTrkAfter(NULL),
290  fHistTPConlyVsCL1Before(NULL),
291  fHistTPConlyVsV0MBefore(NULL),
292  fHistTPConlyVsCL1After(NULL),
293  fHistTPConlyVsV0MAfter(NULL),
294  fHistGlobalVsV0MBefore(NULL),
295  fHistGlobalVsV0MAfter(NULL),
296  fHistRawVsCorrMultFB(NULL),
297  hCentvsTPCmultCuts(NULL),
298  fHV0AEventPlaneVsCent(NULL),
299  fHV0CEventPlaneVsCent(NULL),
300  fHTPCAEventPlaneVsCent(NULL),
301  fHTPCCEventPlaneVsCent(NULL),
302  fHTPCEventPlaneVsCent(NULL),
303  fV0MultChVsRun(NULL),
304  fCentDistBefore(NULL),
305  fCentDistAfter(NULL),
306  fHCorrectV0M(NULL),
307  fHAvgerageQnV0A(NULL),
308  fHAvgerageQnV0C(NULL),
309  fHCentWeightForRun(NULL),
310  fQAEtaPhiAfterNUA(NULL),
311  fQAEtaPhiAfterNUAPion(NULL),
312  fQAEtaPhiAfterNUAKaon(NULL),
313  fQAEtaPhiAfterNUAProton(NULL),
314  fV0AQ2xVsCentRun(NULL),
315  fV0AQ2yVsCentRun(NULL),
316  fV0CQ2xVsCentRun(NULL),
317  fV0CQ2yVsCentRun(NULL),
318  fV0AQ3xVsCentRun(NULL),
319  fV0AQ3yVsCentRun(NULL),
320  fV0CQ3xVsCentRun(NULL),
321  fV0CQ3yVsCentRun(NULL),
322  fTPCAQ2xVsCentRun(NULL),
323  fTPCAQ2yVsCentRun(NULL),
324  fTPCCQ2xVsCentRun(NULL),
325  fTPCCQ2yVsCentRun(NULL),
326  fTPCAQ3xVsCentRun(NULL),
327  fTPCAQ3yVsCentRun(NULL),
328  fTPCCQ3xVsCentRun(NULL),
329  fTPCCQ3yVsCentRun(NULL),
330  fTPCAQ4xVsCentRun(NULL),
331  fTPCAQ4yVsCentRun(NULL),
332  fTPCCQ4xVsCentRun(NULL),
333  fTPCCQ4yVsCentRun(NULL),
334  fTPCFQ2xVsCentRun(NULL),
335  fTPCFQ2yVsCentRun(NULL),
336  fZPASignalPerChVsCent(NULL),
337  fZPCSignalPerChVsCent(NULL),
338  fZNASignalPerChVsCent(NULL),
339  fZNCSignalPerChVsCent(NULL),
340  fCentDistVsVz(NULL),
341  fFilterBit(1),
342  gN(1),
343  gM(1),
344  gPsiN(2),
345  fOldRunNum(111),
346  fEventCount(0),
347  fTPCclustMin(70),
348  fNSigmaCut(2),
349  fMinPtCut(0.2),
350  fMaxPtCut(5.0),
351  fMinEtaCut(-0.8),
352  fMaxEtaCut(0.8),
353  fTrkChi2Min(0.1),
354  fDCAxyMax(2.4),
355  fDCAzMax(3.2),
356  fdEdxMin(10.),
357  fCentralityPercentMin(0),
358  fCentralityPercentMax(90),
359  fPileUpSlopeParm(3.43),
360  fPileUpConstParm(43),
361  fMinVzCut(-10.0),
362  fMaxVzCut(10.0),
363  bApplyMCcorr(kFALSE),
364  bV0MGainCorr(kFALSE),
365  bSkipPileUpCut(kFALSE),
366  bFillNUAHistPID(kFALSE),
367  bUseKinkTracks(kFALSE),
368  sPathOfMCFile("/alien"),
369  sNucleiTP("PbPb"),
370  sCentrEstimator("V0M"),
371  fHistEventCount(NULL)
372 {
373  for(int i=0;i<3;i++){
374  fHistPtwithTPCNsigma[i]=NULL;
375  fHistPtwithTOFmasscut[i]=NULL;
376  fHistPtwithTOFSignal[i]=NULL;
377  fHistTOFnSigmavsPtAfter[i]=NULL;
378  fHistTPCnSigmavsPtAfter[i]=NULL;
380  fHistTPCdEdxvsPtPIDAfter[i]=NULL;
381  }
382  for(int i=0;i<5;i++){
383  fHCorrectNUApos[i] = NULL;
384  fHCorrectNUAneg[i] = NULL;
385  }
386  for(int i=0;i<5;i++){ // for PID NUA
387  fHCorrectNUAposPion[i] = NULL;
388  fHCorrectNUAnegPion[i] = NULL;
389  fHCorrectNUAposKaon[i] = NULL;
390  fHCorrectNUAnegKaon[i] = NULL;
391  fHCorrectNUAposProton[i] = NULL;
392  fHCorrectNUAnegProton[i] = NULL;
393  }
394  //3p vs Centrality
395  for(int i=0;i<2;i++){
396  for(int j=0;j<4;j++){
397  fHist_Corr3p_EP_Norm_PN[i][j] = NULL;
398  fHist_Corr3p_EP_Norm_PP[i][j] = NULL;
399  fHist_Corr3p_EP_Norm_NN[i][j] = NULL;
400  }
401  for(int j=0;j<4;j++) {
402  fHist_Reso2n_EP_Norm_Det[i][j] = NULL;
403  }
404  }
405 
406  //2p vs Centrality:
407 
408  //2p vs Refm:
409 
410  //Differential Charge:
411  for(int i=0;i<2;i++){
412  for(int j=0;j<6;j++){
413  fHist_Corr3p_pTSum_EP_V0A_PN[i][j] = NULL;
414  fHist_Corr3p_pTSum_EP_V0A_PP[i][j] = NULL;
415  fHist_Corr3p_pTSum_EP_V0A_NN[i][j] = NULL;
416  fHist_Corr3p_pTSum_EP_V0C_PN[i][j] = NULL;
417  fHist_Corr3p_pTSum_EP_V0C_PP[i][j] = NULL;
418  fHist_Corr3p_pTSum_EP_V0C_NN[i][j] = NULL;
419 
420  fHist_Corr3p_pTDiff_EP_V0A_PN[i][j] = NULL;
421  fHist_Corr3p_pTDiff_EP_V0A_PP[i][j] = NULL;
422  fHist_Corr3p_pTDiff_EP_V0A_NN[i][j] = NULL;
423  fHist_Corr3p_pTDiff_EP_V0C_PN[i][j] = NULL;
424  fHist_Corr3p_pTDiff_EP_V0C_PP[i][j] = NULL;
425  fHist_Corr3p_pTDiff_EP_V0C_NN[i][j] = NULL;
426 
427  fHist_Corr3p_EtaDiff_EP_V0A_PN[i][j] = NULL;
428  fHist_Corr3p_EtaDiff_EP_V0A_PP[i][j] = NULL;
429  fHist_Corr3p_EtaDiff_EP_V0A_NN[i][j] = NULL;
430  fHist_Corr3p_EtaDiff_EP_V0C_PN[i][j] = NULL;
431  fHist_Corr3p_EtaDiff_EP_V0C_PP[i][j] = NULL;
432  fHist_Corr3p_EtaDiff_EP_V0C_NN[i][j] = NULL;
433  }
434  }
435  for(int i=0;i<4;i++){
436  for(int j=0;j<5;j++){
437  fHist3DEtaPhiVz_Pos_Run[i][j]=NULL;
438  fHist3DEtaPhiVz_Neg_Run[i][j]=NULL;
439  }
440  }
441  for(int i=0;i<10;i++){
442  fFB_Efficiency_Cent[i] = NULL;
443  //fFB_Efficiency_Pion_Cent[i] = NULL;
444  //fFB_Efficiency_Kaon_Cent[i] = NULL;
445  //fFB_Efficiency_Proton_Pos_Cent[i] = NULL;
446  //fFB_Efficiency_Proton_Neg_Cent[i] = NULL;
447  }
448 }
449 
450 //___________________________ destructor ___________________________
452 {
453  //Destructor
454  //if(fPIDResponse) delete fPIDResponse;
455  //if(fMultSelection) delete fMultSelection;
456 
457  /*
458  if(fHCorrectV0M) delete fHCorrectV0M;
459  if(fHAvgerageQnV0A) delete fHAvgerageQnV0A;
460  if(fHAvgerageQnV0C) delete fHAvgerageQnV0C;
461 
462  if(mfileFBHijing->IsOpen()){
463  mfileFBHijing->Close();
464  if(fListFBHijing) delete fListFBHijing;
465  }
466  for(int i=0;i<10;i++){
467  if(fFB_Efficiency_Cent[i])
468  delete fFB_Efficiency_Cent[i];
469  }
470  for(int i=0;i<5;i++){
471  if(fHCorrectNUApos[i]) delete fHCorrectNUApos[i];
472  if(fHCorrectNUAneg[i]) delete fHCorrectNUAneg[i];
473  }
474  for(int i=0;i<5;i++){ // for PID
475  if(fHCorrectNUAposPion[i]) delete fHCorrectNUAposPion[i];
476  if(fHCorrectNUAnegPion[i]) delete fHCorrectNUAnegPion[i];
477  if(fHCorrectNUAposKaon[i]) delete fHCorrectNUAposKaon[i];
478  if(fHCorrectNUAnegKaon[i]) delete fHCorrectNUAnegKaon[i];
479  if(fHCorrectNUAposProton[i]) delete fHCorrectNUAposProton[i];
480  if(fHCorrectNUAnegProton[i]) delete fHCorrectNUAnegProton[i];
481  }
482 */
483 
484  //Delete the clones
485  if(fListFBHijing) delete fListFBHijing;
486  if(fListNUACorr) delete fListNUACorr;
487  if(fListV0MCorr) delete fListV0MCorr;
488 
489  if(fListHist) delete fListHist;
490  if(fAnalysisUtil) delete fAnalysisUtil; // its 'new' !!
491 
492 }//---------------- sanity ------------------------
493 
494 
495 
497 {
498 
499  //std::cout<<"\n..UserCreateOutputObject called.. with isCorr = "<<isCorr<<"\n ....check if succeeded...\n"<<endl;
500 
501  //input hander
502  AliAnalysisManager *mgr=AliAnalysisManager::GetAnalysisManager();
503  AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(mgr->GetInputEventHandler());
504  if (!inputHandler) { printf("\n ...Input handler missing!!!...\n"); return; }
505 
506  //PileUp Multi-Vertex
507  fAnalysisUtil = new AliAnalysisUtils();
508  fAnalysisUtil->SetUseMVPlpSelection(kTRUE);
509  fAnalysisUtil->SetUseOutOfBunchPileUp(kTRUE);
510 
511  //pid response object
512  //fPIDResponse=inputHandler->GetPIDResponse();
513 
514  fListHist = new TList();
515  fListHist->SetOwner(kTRUE);
516 
517 
519 
520  if(!gGrid){
521  TGrid::Connect("alien://");
522  }
523 
525 
526 
527 
528  hCentvsTPCmultCuts = new TH2F("hCentvsTPCmultCuts","TPCmult Low,high",100,0,100,5,0,5);
529  hCentvsTPCmultCuts->GetYaxis()->SetBinLabel(1,"mean");
530  hCentvsTPCmultCuts->GetYaxis()->SetBinLabel(2,"sigma");
532 
534 
535  fHistTOFMatchCount = new TH2F("fHistTOFMatchCount","TofMatchFlag vs mismatch Prob",10,0,10,200,-5,5);
537 
538  fHistEtaPtBefore = new TH2F("fHistEtaPtBefore","Eta vs pT",100,-1.25,1.25,100,0,10);
540 
541  fHistEtaPtAfter = new TH2F("fHistEtaPtAfter","Eta vs pT",100,-1.25,1.25,100,0,10);
543 
544 
545  Int_t gMaxTPCFB1mult = 0;
546  Int_t gMaxGlobalmult = 0;
547  Int_t gMaxTPCcorrmult = 0;
548  Int_t gMaxESDtracks = 0;
549  Int_t nBinRefMult = 200;
550  Int_t nRefMultMax = 4000;
551 
552 
553 
554  if(sNucleiTP=="pp"||sNucleiTP=="PP"){
555  gMaxGlobalmult = 200;
556  gMaxTPCFB1mult = 200;
557  gMaxTPCcorrmult = 500;
558  gMaxESDtracks = 1000;
559  nBinRefMult = 100; //change binning for pp
560  nRefMultMax = 500;
561  //fSkipOutlierCut = 1;
562  }
563  else if(sNucleiTP=="pPb"||sNucleiTP=="Pbp"||sNucleiTP=="PbP"||sNucleiTP=="PPb"){
564  gMaxGlobalmult = 400;
565  gMaxTPCFB1mult = 400;
566  gMaxTPCcorrmult = 500;
567  gMaxESDtracks = 2000;
568  nBinRefMult = 100; //change binning for pPb
569  nRefMultMax = 1000;
570  //fSkipOutlierCut = 1;
571  }
572  else{
573  gMaxGlobalmult = 4000;
574  gMaxTPCFB1mult = 4000;
575  gMaxTPCcorrmult = 5000;
576  gMaxESDtracks = 20000;
577  //fSkipOutlierCut = 0;
578  }
579 
580 
581 
582  //if(bSkipPileUpCut) { fSkipOutlierCut = 1;}
583 
584 
585  fHistTPCVsESDTrkBefore = new TH2F("fHistTPCVsESDTrkBefore","Before; TPC1; ESD trk",500,0,gMaxTPCcorrmult,200,0,gMaxESDtracks);
587  fHistTPCVsESDTrkAfter = new TH2F("fHistTPCVsESDTrkAfter"," After; TPC1; ESD trk",500,0,gMaxTPCcorrmult,200,0,gMaxESDtracks);
589 
590  fHistTPCvsGlobalMultBefore = new TH2F("fHistTPCvsGlobalMultBefore","Before; Global; TPC(fb1) ",200,0,gMaxGlobalmult,200,0,gMaxTPCFB1mult);
592  fHistTPCvsGlobalMultAfter = new TH2F("fHistTPCvsGlobalMultAfter"," After; Global; TPC(fb1) ",200,0,gMaxGlobalmult,200,0,gMaxTPCFB1mult);
594 
595  fHistGlobalVsV0MBefore = new TH2F("fHistGlobalVsV0MBefore","Before;Cent(V0M);Global",100,0,100,500,0,gMaxGlobalmult);
597  fHistGlobalVsV0MAfter = new TH2F("fHistGlobalVsV0MAfter"," After; Cent(V0M);Global",100,0,100,500,0,gMaxGlobalmult);
599 
600  fHistTPConlyVsCL1Before = new TH2F("fHistTPConlyVsCL1Before","Before;Cent(CL1); TPC(FB1)",100,0,100,250,0,gMaxTPCcorrmult);
602  fHistTPConlyVsCL1After = new TH2F("fHistTPConlyVsCL1After","After; Cent(CL1); TPC(FB1) ",100,0,100,250,0,gMaxTPCcorrmult);
604 
605  fHistTPConlyVsV0MBefore = new TH2F("fHistTPConlyVsV0MBefore","Before;Cent(V0M); TPC(FB1)",100,0,100,250,0,gMaxTPCcorrmult);
607  fHistTPConlyVsV0MAfter = new TH2F("fHistTPConlyVsV0MAfter","After; Cent(V0M); TPC(FB1) ",100,0,100,250,0,gMaxTPCcorrmult);
609 
610  fHistRawVsCorrMultFB = new TH2F("fHistRawVsCorrMultFB",Form("FB%d;Mult_{raw};Mult_{corr}",fFilterBit),gMaxTPCFB1mult,0,gMaxTPCFB1mult,gMaxTPCcorrmult,0,gMaxTPCcorrmult);
612 
613 
614 
615 
616 
617 
618  // Turn off PID QA histograms
619  //---------------- PID QA Histograms ---------------------
620  /*
621  fHistTPCdEdxvsPBefore = new TH2F("fHistTPCdEdxvsPBefore","Before; p (GeV/c); dEdx (arb)",200,-5,5,200,0,250);
622  fListHist->Add(fHistTPCdEdxvsPBefore);
623  fHistTPCdEdxvsPAfter = new TH2F("fHistTPCdEdxvsPAfter"," After; p (GeV/c); dEdx (arb)",200,-5,5, 200,0,250);
624  fListHist->Add(fHistTPCdEdxvsPAfter);
625 
626  fHistTOFBetavsPBefore = new TH2F("fHistTOFBetavsPBefore","Before; p (GeV/c); beta ",200,-5,5,100,0.0,1.2);
627  fListHist->Add(fHistTOFBetavsPBefore);
628  fHistTOFBetavsPAfter = new TH2F("fHistTOFBetavsPAfter"," After; p (GeV/c); beta ",200,-5,5,100,0.0,1.2);
629  fListHist->Add(fHistTOFBetavsPAfter);
630 
631  fHistTOFMassvsPtBefore = new TH2F("fHistTOFMassvsPtBefore","Before; p_{T}(GeV/c); m^{2}(GeV^{2}/c^{4})",200,-5,5,500,-0.5,4.5);
632  fListHist->Add(fHistTOFMassvsPtBefore);
633 
634  //const char *gSpecies[4] = {"Pion","Kaon","proton","Charge"};
635 
636  for(int i=0;i<3;i++){
637  fHistPtwithTPCNsigma[i] = new TH1F(Form("fHistPtwithTPCNsigma_%d",i), Form("%d;p_{T}(GeV/c))",i),200,-5,5); //i
638  fListHist->Add(fHistPtwithTPCNsigma[i]);
639  fHistPtwithTOFmasscut[i] = new TH1F(Form("fHistPtwithTOFmasscut_%d",i),Form("%d;p_{T}(GeV/c))",i),200,-5,5);
640  fListHist->Add(fHistPtwithTOFmasscut[i]);
641  fHistPtwithTOFSignal[i] = new TH1F(Form("fHistPtwithTOFSignal_%d", i),Form("%d;p_{T}(GeV/c))",i),200,-5,5);
642  fListHist->Add(fHistPtwithTOFSignal[i]);
643 
644  fHistTOFnSigmavsPtAfter[i] = new TH2F(Form("fHistTOFnSigmavsPtAfter_%d",i),Form("%d;p_{T}(GeV/c);n#sigma_{TOF}",i),200,-5,5,400,-10.0,10.0);
645  fListHist->Add(fHistTOFnSigmavsPtAfter[i]);
646 
647  fHistTPCnSigmavsPtAfter[i] = new TH2F(Form("fHistTPCnSigmavsPtAfter_%d",i),Form("%d;p_{T}(GeV/c);n#sigma_{TPC}",i),200,-5,5,400,-10.0,10.0);
648  fListHist->Add(fHistTPCnSigmavsPtAfter[i]);
649 
650  fHistTPCTOFnSigmavsPtAfter[i] = new TH3F(Form("fHistTPCTOFnSigmavsPtAfter_%d",i),Form("%d; p_{T}(GeV/c); n#sigma_{TPC}; n#sigma_{TOF}",i),100,0,5,400,-10,10,400,-10,10);
651  fListHist->Add(fHistTPCTOFnSigmavsPtAfter[i]);
652 
653  fHistTPCdEdxvsPtPIDAfter[i] = new TH2F(Form("fHistTPCdEdxvsPtAfter_%d",i),"AfterCut; p_{T} (GeV/c); dEdx (arb)",400,0,10,200,0,250);
654  fListHist->Add(fHistTPCdEdxvsPtPIDAfter[i]);
655  }// PID histograms done
656  */
657 
658 
659 
660 
661  fCentDistBefore = new TH1F("fCentDistBefore","no Cut; Cent (%); Events ",100,0,100);
663 
664  fCentDistAfter = new TH1F("fCentDistAfter","with Cut; Cent (%); Events ",100,0,100);
666 
667 
668 
669 
670 
671 
672 
673 
674 
675 
676 
677 
678 
679  Double_t centRange[11] = {0,5,10,20,30,40,50,60,70,80,90};
680  //const char *gDetForEP[4] = {"V0A","V0C","TPC-A","TPC-C"};
681  // 10,centRange
682  //------------------- 3p correlator vs Centrality (EP method) ------------------
683  for(int i=0;i<2;i++){
684  //Charged:
685  for(int j=0;j<4;j++){
686  //Detector: 0 = V0A, 1 = V0C, 3 = TPCA, 4 = TPCC
687  fHist_Corr3p_EP_Norm_PN[i][j] = new TProfile(Form("fHist_Corr3p_EP_Norm_PosNeg_Mag%d_Det%d",i,j+1),Form("US, #Psi_{2} %d",j),10,centRange,"");
688  fHist_Corr3p_EP_Norm_PN[i][j]->Sumw2();
690  fHist_Corr3p_EP_Norm_PP[i][j] = new TProfile(Form("fHist_Corr3p_EP_Norm_PosPos_Mag%d_Det%d",i,j+1),Form("P-P, #Psi_{2} %d",j),10,centRange,"");
691  fHist_Corr3p_EP_Norm_PP[i][j]->Sumw2();
693  fHist_Corr3p_EP_Norm_NN[i][j] = new TProfile(Form("fHist_Corr3p_EP_Norm_NegNeg_Mag%d_Det%d",i,j+1),Form("N-N, #Psi_{2}, %d",j),10,centRange,"");
694  fHist_Corr3p_EP_Norm_NN[i][j]->Sumw2();
696  }
697  //EP Resolution:
698  for(int j=0;j<4;j++){
699  //Det: 0 = v0c-v0a, 1 = v0a-TPC, 2 = v0c-TPC, 3 =TPC-A TPC-C
700  fHist_Reso2n_EP_Norm_Det[i][j] = new TProfile(Form("fHist_Reso2n_EP_Norm_Mag%d_DetComb%d",i,j+1),"Event plane Resolution",10,centRange,"");
701  fHist_Reso2n_EP_Norm_Det[i][j]->Sumw2();
703  }
704  }//magfield loop
705 
706 
707 
708 
709  Char_t name[100];
710  Char_t title[100];
711 
712  Double_t pTRange[21] = {0.0,0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,1.8,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0};
713 //Double_t pTRange[11] = {0.0,0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6,1.8,2.0};
714 
715  //Charge:
716  for(Int_t i=0;i<2;i++){
717  for(Int_t j=0;j<6;j++){
718  sprintf(name,"fHist_Corr3p_pTSum_EP_V0A_PN_Mag%d_Cent%d",i,j);
719  sprintf(title,"PN 3p vs (pT1+pT2)/2, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
720  fHist_Corr3p_pTSum_EP_V0A_PN[i][j] = new TProfile(name,title,20,pTRange,"");
721  fHist_Corr3p_pTSum_EP_V0A_PN[i][j]->Sumw2();
723 
724  sprintf(name,"fHist_Corr3p_pTSum_EP_V0A_PP_Mag%d_Cent%d",i,j);
725  sprintf(title,"PP 3p vs (pT1+pT2)/2, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
726  fHist_Corr3p_pTSum_EP_V0A_PP[i][j] = new TProfile(name,title,20,pTRange,"");
727  fHist_Corr3p_pTSum_EP_V0A_PP[i][j]->Sumw2();
729 
730  sprintf(name,"fHist_Corr3p_pTSum_EP_V0A_NN_Mag%d_Cent%d",i,j);
731  sprintf(title,"NN 3p vs (pT1+pT2)/2, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
732  fHist_Corr3p_pTSum_EP_V0A_NN[i][j] = new TProfile(name,title,20,pTRange,"");
733  fHist_Corr3p_pTSum_EP_V0A_NN[i][j]->Sumw2();
735  //-----v0c----
736  sprintf(name,"fHist_Corr3p_pTSum_EP_V0C_PN_Mag%d_Cent%d",i,j);
737  sprintf(title,"PN 3p vs (pT1+pT2)/2, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
738  fHist_Corr3p_pTSum_EP_V0C_PN[i][j] = new TProfile(name,title,20,pTRange,"");
739  fHist_Corr3p_pTSum_EP_V0C_PN[i][j]->Sumw2();
741 
742  sprintf(name,"fHist_Corr3p_pTSum_EP_V0C_PP_Mag%d_Cent%d",i,j);
743  sprintf(title,"PP 3p vs (pT1+pT2)/2, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
744  fHist_Corr3p_pTSum_EP_V0C_PP[i][j] = new TProfile(name,title,20,pTRange,"");
745  fHist_Corr3p_pTSum_EP_V0C_PP[i][j]->Sumw2();
747 
748  sprintf(name,"fHist_Corr3p_pTSum_EP_V0C_NN_Mag%d_Cent%d",i,j);
749  sprintf(title,"NN 3p vs (pT1+pT2)/2, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
750  fHist_Corr3p_pTSum_EP_V0C_NN[i][j] = new TProfile(name,title,20,pTRange,"");
751  fHist_Corr3p_pTSum_EP_V0C_NN[i][j]->Sumw2();
753  }
754  }
755 
756  for(Int_t i=0;i<2;i++){
757  for(Int_t j=0;j<6;j++){
758  sprintf(name,"fHist_Corr3p_pTDiff_EP_V0A_PN_Mag%d_Cent%d",i,j);
759  sprintf(title,"PN 3p vs |pT1-pT2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
760  fHist_Corr3p_pTDiff_EP_V0A_PN[i][j] = new TProfile(name,title,20,pTRange,"");
761  fHist_Corr3p_pTDiff_EP_V0A_PN[i][j]->Sumw2();
763 
764  sprintf(name,"fHist_Corr3p_pTDiff_EP_V0A_PP_Mag%d_Cent%d",i,j);
765  sprintf(title,"PP 3p vs |pT1-pT2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
766  fHist_Corr3p_pTDiff_EP_V0A_PP[i][j] = new TProfile(name,title,20,pTRange,"");
767  fHist_Corr3p_pTDiff_EP_V0A_PP[i][j]->Sumw2();
769 
770  sprintf(name,"fHist_Corr3p_pTDiff_EP_V0A_NN_Mag%d_Cent%d",i,j);
771  sprintf(title,"NN 3p vs |pT1-pT2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
772  fHist_Corr3p_pTDiff_EP_V0A_NN[i][j] = new TProfile(name,title,20,pTRange,"");
773  fHist_Corr3p_pTDiff_EP_V0A_NN[i][j]->Sumw2();
775  //-----v0c----
776  sprintf(name,"fHist_Corr3p_pTDiff_EP_V0C_PN_Mag%d_Cent%d",i,j);
777  sprintf(title,"PN 3p vs |pT1-pT2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
778  fHist_Corr3p_pTDiff_EP_V0C_PN[i][j] = new TProfile(name,title,20,pTRange,"");
779  fHist_Corr3p_pTDiff_EP_V0C_PN[i][j]->Sumw2();
781 
782  sprintf(name,"fHist_Corr3p_pTDiff_EP_V0C_PP_Mag%d_Cent%d",i,j);
783  sprintf(title,"PP 3p vs |pT1-pT2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
784  fHist_Corr3p_pTDiff_EP_V0C_PP[i][j] = new TProfile(name,title,20,pTRange,"");
785  fHist_Corr3p_pTDiff_EP_V0C_PP[i][j]->Sumw2();
787 
788  sprintf(name,"fHist_Corr3p_pTDiff_EP_V0C_NN_Mag%d_Cent%d",i,j);
789  sprintf(title,"NN 3p vs |pT1-pT2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
790  fHist_Corr3p_pTDiff_EP_V0C_NN[i][j] = new TProfile(name,title,20,pTRange,"");
791  fHist_Corr3p_pTDiff_EP_V0C_NN[i][j]->Sumw2();
793  }
794  }
795 
796  //Double_t EtaRange[9] = {0.0,0.2,0.4,0.6,0.8,1.0,1.2,1.4,1.6}; //Use this after tests done
797  //Now Eta binning: 16,0,1.6 for test
798 
799  for(Int_t i=0;i<2;i++){
800  for(Int_t j=0;j<6;j++){
801  sprintf(name,"fHist_Corr3p_EtaDiff_EP_V0A_PN_Mag%d_Cent%d",i,j);
802  sprintf(title,"PN 3p vs |Eta1-Eta2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
803  fHist_Corr3p_EtaDiff_EP_V0A_PN[i][j] = new TProfile(name,title,8,0,1.6,"");
804  fHist_Corr3p_EtaDiff_EP_V0A_PN[i][j]->Sumw2();
806 
807  sprintf(name,"fHist_Corr3p_EtaDiff_EP_V0A_PP_Mag%d_Cent%d",i,j);
808  sprintf(title,"PP 3p vs |Eta1-Eta2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
809  fHist_Corr3p_EtaDiff_EP_V0A_PP[i][j] = new TProfile(name,title,8,0,1.6,"");
810  fHist_Corr3p_EtaDiff_EP_V0A_PP[i][j]->Sumw2();
812 
813  sprintf(name,"fHist_Corr3p_EtaDiff_EP_V0A_NN_Mag%d_Cent%d",i,j);
814  sprintf(title,"NN 3p vs |Eta1-Eta2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
815  fHist_Corr3p_EtaDiff_EP_V0A_NN[i][j] = new TProfile(name,title,8,0,1.6,"");
816  fHist_Corr3p_EtaDiff_EP_V0A_NN[i][j]->Sumw2();
818  //-----v0c----
819  sprintf(name,"fHist_Corr3p_EtaDiff_EP_V0C_PN_Mag%d_Cent%d",i,j);
820  sprintf(title,"PN 3p vs |Eta1-Eta2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
821  fHist_Corr3p_EtaDiff_EP_V0C_PN[i][j] = new TProfile(name,title,8,0,1.6,"");
822  fHist_Corr3p_EtaDiff_EP_V0C_PN[i][j]->Sumw2();
824 
825  sprintf(name,"fHist_Corr3p_EtaDiff_EP_V0C_PP_Mag%d_Cent%d",i,j);
826  sprintf(title,"PP 3p vs |Eta1-Eta2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
827  fHist_Corr3p_EtaDiff_EP_V0C_PP[i][j] = new TProfile(name,title,8,0,1.6,"");
828  fHist_Corr3p_EtaDiff_EP_V0C_PP[i][j]->Sumw2();
830 
831  sprintf(name,"fHist_Corr3p_EtaDiff_EP_V0C_NN_Mag%d_Cent%d",i,j);
832  sprintf(title,"NN 3p vs |Eta1-Eta2|, Cent %2.0f-%2.0f",centRange[i],centRange[i+1]);
833  fHist_Corr3p_EtaDiff_EP_V0C_NN[i][j] = new TProfile(name,title,8,0,1.6,"");
834  fHist_Corr3p_EtaDiff_EP_V0C_NN[i][j]->Sumw2();
836  }
837  }
838 
839 
840 
841  //---- to store NUA and calib histograms -----
842  TList *fListNUACalib = new TList();
843  fListNUACalib->SetName("fListNUACalib");
844  fListNUACalib->SetOwner(kTRUE);
845 
846 
847  Double_t truncPi = 3.1416;
848 
849  //-------------- Event plane distributions --------------
850  fHV0AEventPlaneVsCent = new TH2F("fHV0AEventPlaneVsCent",Form("Psi %d from V0A", gPsiN), 10,centRange,50,-0.0,truncPi);
851  fListNUACalib->Add(fHV0AEventPlaneVsCent);
852 
853  fHV0CEventPlaneVsCent = new TH2F("fHV0CEventPlaneVsCent",Form("Psi %d from V0C", gPsiN), 10,centRange,50,-0.0,truncPi);
854  fListNUACalib->Add(fHV0CEventPlaneVsCent);
855 
856  fHTPCAEventPlaneVsCent = new TH2F("fHTPCAEventPlaneVsCent",Form("Psi %d from pos eta",gPsiN),10,centRange,50,-0.0,truncPi);
857  fListNUACalib->Add(fHTPCAEventPlaneVsCent);
858 
859  fHTPCCEventPlaneVsCent = new TH2F("fHTPCCEventPlaneVsCent",Form("Psi %d from neg eta",gPsiN),10,centRange,50,-0.0,truncPi);
860  fListNUACalib->Add(fHTPCCEventPlaneVsCent);
861 
862  fHTPCEventPlaneVsCent = new TH2F("fHTPCEventPlaneVsCent",Form("Psi %d from Full TPC",gPsiN),10,centRange,50,-0.0,truncPi);
863  fListNUACalib->Add(fHTPCEventPlaneVsCent);
864 
865 
866 
867  //-------------------- QA: how does the eta-phi looks after NUA correction ----------------
868  fQAEtaPhiAfterNUA = new TH2F("fQAEtaPhiAfterNUA","eta vs phi with NUA corr",50,0,6.283185,16,-0.8,0.8);
869  fListNUACalib->Add(fQAEtaPhiAfterNUA);
870 
871  fQAEtaPhiAfterNUAPion = new TH2F("fQAEtaPhiAfterNUAPion","Pion eta vs phi with NUA corr",50,0,6.283185,16,-0.8,0.8);
872  fListNUACalib->Add(fQAEtaPhiAfterNUAPion);
873 
874  fQAEtaPhiAfterNUAKaon = new TH2F("fQAEtaPhiAfterNUAKaon","Kaon eta vs phi with NUA corr",50,0,6.283185,16,-0.8,0.8);
875  fListNUACalib->Add(fQAEtaPhiAfterNUAKaon);
876 
877  fQAEtaPhiAfterNUAProton = new TH2F("fQAEtaPhiAfterNUAProton","Proton eta vs phi with NUA corr",50,0,6.283185,16,-0.8,0.8);
878  fListNUACalib->Add(fQAEtaPhiAfterNUAProton);
879 
880 
881 
882  //----------------- ZPA,ZNA Calibration hist: ---------------------
883  fZPASignalPerChVsCent = new TH2F("fZPASignalPerChVsCent","En-ZPA, Cent, Vz, ch",91,0,91,5,0,5);
884  fZPCSignalPerChVsCent = new TH2F("fZPCSignalPerChVsCent","En-ZPC, Cent, Vz, ch",91,0,91,5,0,5);
885  fZNASignalPerChVsCent = new TH2F("fZNASignalPerChVsCent","En-ZNA, Cent, Vz, ch",91,0,91,5,0,5);
886  fZNCSignalPerChVsCent = new TH2F("fZNCSignalPerChVsCent","En-ZNC, Cent, Vz, ch",91,0,91,5,0,5);
887 
888  fListNUACalib->Add(fZPASignalPerChVsCent);
889  fListNUACalib->Add(fZPCSignalPerChVsCent);
890  fListNUACalib->Add(fZNASignalPerChVsCent);
891  fListNUACalib->Add(fZNCSignalPerChVsCent);
892 
893  fCentDistVsVz = new TProfile("fCentDistVsVz","<Cent> vs Vz.",80,-10,10);
894 
895  fListNUACalib->Add(fCentDistVsVz);
896 
897 
898 
899 
900  //----------------- V0 Calibration hist: ---------------------
901  fV0MultChVsRun = new TH2F("fV0MultChVsRun","1-32 V0C, 33-64 V0A",64,0,64,90,0,90);
902  fListNUACalib->Add(fV0MultChVsRun);
903 
904  //const char *sCorrect[2]={"wo Corr","w Corr"};
905  Int_t isCorr = 0;
906 
907  if(fListV0MCorr){
908  isCorr = 1;
909  }
910 
911 
912  fV0AQ2xVsCentRun = new TProfile("fV0ACos2nVsCentRun",Form("<Cos2> vs cent (%d)",isCorr),90,0,90,""); //sCorrect[isCorr]
913  fListNUACalib->Add(fV0AQ2xVsCentRun);
914  fV0AQ2yVsCentRun = new TProfile("fV0ASin2nVsCentRun",Form("<Sin2> vs cent (%d)",isCorr),90,0,90,"");
915  fListNUACalib->Add(fV0AQ2yVsCentRun);
916  fV0CQ2xVsCentRun = new TProfile("fV0CCos2nVsCentRun",Form("<Cos2> vs cent (%d)",isCorr),90,0,90,"");
917  fListNUACalib->Add(fV0CQ2xVsCentRun);
918  fV0CQ2yVsCentRun = new TProfile("fV0CSin2nVsCentRun",Form("<Sin2> vs cent (%d)",isCorr),90,0,90,"");
919  fListNUACalib->Add(fV0CQ2yVsCentRun);
920 
921  fV0AQ3xVsCentRun = new TProfile("fV0ACos3nVsCentRun",Form("<Cos3> vs cent (%d)",isCorr),90,0,90,"");
922  fListNUACalib->Add(fV0AQ3xVsCentRun);
923  fV0AQ3yVsCentRun = new TProfile("fV0ASin3nVsCentRun",Form("<Sin3> vs cent (%d)",isCorr),90,0,90,"");
924  fListNUACalib->Add(fV0AQ3yVsCentRun);
925  fV0CQ3xVsCentRun = new TProfile("fV0CCos3nVsCentRun",Form("<Cos3> vs cent (%d)",isCorr),90,0,90,"");
926  fListNUACalib->Add(fV0CQ3xVsCentRun);
927  fV0CQ3yVsCentRun = new TProfile("fV0CSin3nVsCentRun",Form("<Sin3> vs cent (%d)",isCorr),90,0,90,"");
928  fListNUACalib->Add(fV0CQ3yVsCentRun);
929 
930  isCorr = 1;
931  if(fListNUACorr){
932  cout<<"\n =========> NUA file found for NUA correction <========== \n";
933  isCorr = 0;
934  }
935  //------------------- TPC Qvector Recentering Histograms --------------
936  fTPCAQ2xVsCentRun = new TProfile("fTPCACos2nVsCentRun",Form("<Cos2> vs cent (%d)",isCorr),90,0,90,"");
937  fListNUACalib->Add(fTPCAQ2xVsCentRun);
938  fTPCAQ2yVsCentRun = new TProfile("fTPCASin2nVsCentRun",Form("<Sin2> vs cent (%d)",isCorr),90,0,90,"");
939  fListNUACalib->Add(fTPCAQ2yVsCentRun);
940  fTPCCQ2xVsCentRun = new TProfile("fTPCCCos2nVsCentRun",Form("<Cos2> vs cent (%d)",isCorr),90,0,90,"");
941  fListNUACalib->Add(fTPCCQ2xVsCentRun);
942  fTPCCQ2yVsCentRun = new TProfile("fTPCCSin2nVsCentRun",Form("<Sin2> vs cent (%d)",isCorr),90,0,90,"");
943  fListNUACalib->Add(fTPCCQ2yVsCentRun);
944 
945  fTPCAQ3xVsCentRun = new TProfile("fTPCACos3nVsCentRun",Form("<Cos3> vs cent (%d)",isCorr),90,0,90,"");
946  fListNUACalib->Add(fTPCAQ3xVsCentRun);
947  fTPCAQ3yVsCentRun = new TProfile("fTPCASin3nVsCentRun",Form("<Sin3> vs cent (%d)",isCorr),90,0,90,"");
948  fListNUACalib->Add(fTPCAQ3yVsCentRun);
949  fTPCCQ3xVsCentRun = new TProfile("fTPCCCos3nVsCentRun",Form("<Cos3> vs cent (%d)",isCorr),90,0,90,"");
950  fListNUACalib->Add(fTPCCQ3xVsCentRun);
951  fTPCCQ3yVsCentRun = new TProfile("fTPCCSin3nVsCentRun",Form("<Sin3> vs cent (%d)",isCorr),90,0,90,"");
952  fListNUACalib->Add(fTPCCQ3yVsCentRun);
953 
954  fTPCAQ4xVsCentRun = new TProfile("fTPCACos4nVsCentRun",Form("<Cos4> vs cent (%d)",isCorr),90,0,90,"");
955  fListNUACalib->Add(fTPCAQ4xVsCentRun);
956  fTPCAQ4yVsCentRun = new TProfile("fTPCASin4nVsCentRun",Form("<Sin4> vs cent (%d)",isCorr),90,0,90,"");
957  fListNUACalib->Add(fTPCAQ4yVsCentRun);
958  fTPCCQ4xVsCentRun = new TProfile("fTPCCCos4nVsCentRun",Form("<Cos4> vs cent (%d)",isCorr),90,0,90,"");
959  fListNUACalib->Add(fTPCCQ4xVsCentRun);
960  fTPCCQ4yVsCentRun = new TProfile("fTPCCSin4nVsCentRun",Form("<Sin4> vs cent (%d)",isCorr),90,0,90,"");
961  fListNUACalib->Add(fTPCCQ4yVsCentRun);
962 
963 
964  fTPCFQ2xVsCentRun = new TProfile("fTPCFQ2xVsCentRun",Form("<Cos2> vs cent (%d)",isCorr),90,0,90,"");
965  fListNUACalib->Add(fTPCFQ2xVsCentRun);
966  fTPCFQ2yVsCentRun = new TProfile("fTPCFQ2yVsCentRun",Form("<Sin2> vs cent (%d)",isCorr),90,0,90,"");
967  fListNUACalib->Add(fTPCFQ2yVsCentRun);
968 
969  //-------------------------------------------------------------------------------
970 
971 
972 
973 
974  //-------------------------- Define NUA Hist for PID -----------------------------
975  Int_t gCentForNUA[6] = {0,5,10,20,40,90};
976  //Char_t name[100];
977  //Char_t title[100];
978 
979  for(int i=0;i<4;i++){
980  for(int j=0;j<5;j++){
981  sprintf(name,"fHistEtaPhiVz_%d_Pos_Cent%d_Run%d",i,j,1); //gSpecies[i]
982  sprintf(title,"eta,phi,Vz %dPos, Cent%d-%d, FB %d",i,gCentForNUA[j],gCentForNUA[j+1],fFilterBit);
983  fHist3DEtaPhiVz_Pos_Run[i][j] = new TH3F(name,title,10,-10,10,50,0,6.283185,16,-0.8,0.8);
984  fListNUACalib->Add(fHist3DEtaPhiVz_Pos_Run[i][j]);
985 
986  sprintf(name,"fHistEtaPhiVz_%d_Neg_Cent%d_Run%d",i,j,1); //gSpecies[i]
987  sprintf(title,"eta,phi,Vz %dNeg, Cent%d-%d, FB %d",i,gCentForNUA[j],gCentForNUA[j+1],fFilterBit);
988  fHist3DEtaPhiVz_Neg_Run[i][j] = new TH3F(name,title,10,-10,10,50,0,6.283185,16,-0.8,0.8);
989  fListNUACalib->Add(fHist3DEtaPhiVz_Neg_Run[i][j]);
990  }
991  }
992  //---------------------------------------------------------------------------------
993 
994 
995  fListHist->Add(fListNUACalib);
996 
997  PostData(1,fListHist);
998  cout<<"\n.........UserCreateOutputObject called.........\n fFilterBit = "<<fFilterBit<<" CentMax = "<<fCentralityPercentMax;
999  cout<<" PU C = "<<fPileUpConstParm<<" gN = "<<gN<<" gM = "<<gM<<" PsiN = "<<gPsiN<<"\n\n"<<endl;
1000  //cout<<" *********** checking my commit...!! ***************** \n\n ";
1001 }
1002 
1003 
1004 
1005 
1006 
1007 
1008 
1009 
1010 
1011 
1012 
1013 
1014 
1015 
1016 
1017 
1018 
1019 //______________________________________________________________________
1021  //debug only
1022  //cout<<"\n Info:UserExec() called ..!!!\n";
1023  //watch.Start(kTRUE);
1024  //if(fEventCount==501) return;
1025 
1026 
1027  Float_t stepCount = 0.5;
1028 
1029  fHistEventCount->Fill(stepCount); //1
1030  stepCount++;
1031 
1032  fAOD = dynamic_cast<AliAODEvent*>(InputEvent());
1033  fESD = dynamic_cast<AliESDEvent*>(InputEvent());
1034 
1035  if(!(fESD || fAOD)){ printf("ERROR: fESD & fAOD not available\n"); return; }
1036 
1037  fVevent = dynamic_cast<AliVEvent*>(InputEvent());
1038 
1039  if (!fVevent) { printf("ERROR: fVevent not available\n"); return; }
1040 
1041  fHistEventCount->Fill(stepCount); //2
1042  stepCount++;
1043 
1044 
1045 
1046 
1047 
1048 
1049  //--------- Check if I have PID response object --------
1050  /*
1051  if(!fPIDResponse){
1052  AliAnalysisManager *mgr=AliAnalysisManager::GetAnalysisManager();
1053  AliInputEventHandler *inputHandler=dynamic_cast<AliInputEventHandler*>(mgr->GetInputEventHandler());
1054  if(inputHandler) fPIDResponse=inputHandler->GetPIDResponse();
1055  if(!fPIDResponse){
1056  printf("\n\n...... PIDResponse object not found..... \n\n");
1057  return;
1058  }
1059  }
1060  */
1061 
1062 
1063 
1064 
1065  //-------------- Vtx cuts ---------------
1066  const AliVVertex *pVtx = fVevent->GetPrimaryVertex();
1067 
1068  Double_t pVtxZ = -999;
1069  pVtxZ = pVtx->GetZ();
1070 
1071 //if(TMath::Abs(pVtxZ)>10.) return;
1072  //User defined cut:
1073  if(pVtxZ<fMinVzCut || pVtxZ>fMaxVzCut ) return;
1074 
1075 
1076  fHistEventCount->Fill(stepCount); //3
1077  stepCount++;
1078 
1079 
1080 
1081  Float_t centrality = -99.0;
1082  Float_t centrV0M = -99.0;
1083  Float_t centrCL1 = -99.0;
1084 
1085 
1086 
1087  //---------- Centrality Estimators -------------
1088  AliCentrality* Alicentr = ((AliVAODHeader*)fAOD->GetHeader())->GetCentralityP(); // For Run1, 2010 data
1089 
1090  fMultSelection = (AliMultSelection*) InputEvent()->FindListObject("MultSelection"); // Must never comment this
1091  if(!fMultSelection) { printf("\n\n **WARNING** \n::UserExec() AliMultSelection object not found.\n\n"); exit(1); }
1092 
1093 
1094  if(sNucleiTP=="PbPb2010" || sNucleiTP=="2010") {
1095  if(Alicentr){
1096 
1097  centrV0M = Alicentr->GetCentralityPercentile("V0M");
1098  centrCL1 = Alicentr->GetCentralityPercentile("CL1");
1099 
1100  if(sCentrEstimator=="V0M" || sCentrEstimator=="V0"){
1101  centrality = centrV0M;
1102  }
1103  else if(sCentrEstimator=="CL1"){
1104  centrality = centrCL1;
1105  }
1106  else if(sCentrEstimator=="V0C"){
1107  centrality = Alicentr->GetCentralityPercentile("V0C");
1108  }
1109  else if(sCentrEstimator=="V0A"){
1110  centrality = Alicentr->GetCentralityPercentile("V0A");
1111  }
1112  else if(sCentrEstimator=="TRK"){
1113  centrality = Alicentr->GetCentralityPercentile("TRK");
1114  }
1115  }
1116  }
1117  else{ // fall back to MultSelection if other than 2010 data
1118 
1119  centrV0M = fMultSelection->GetMultiplicityPercentile("V0M");
1120  centrCL1 = fMultSelection->GetMultiplicityPercentile("CL1");
1121 
1122  if(sCentrEstimator=="V0M" || sCentrEstimator=="V0"){
1123  centrality = centrV0M;
1124  }
1125  else if(sCentrEstimator=="CL1"){
1126  centrality = centrCL1;
1127  }
1128  else if(sCentrEstimator=="V0C"){
1129  centrality = fMultSelection->GetMultiplicityPercentile("V0C");
1130  }
1131  else if(sCentrEstimator=="V0A"){
1132  centrality = fMultSelection->GetMultiplicityPercentile("V0A");
1133  }
1134  else if(sCentrEstimator=="TRK"){
1135  centrality = fMultSelection->GetMultiplicityPercentile("TRK");
1136  }
1137  }
1138 
1139  /*
1140  //--------- Centrality Estimators -------------
1141  centrV0M = fMultSelection->GetMultiplicityPercentile("V0M");
1142  centrCL1 = fMultSelection->GetMultiplicityPercentile("CL1");
1143 
1144  if(sCentrEstimator=="V0M" || sCentrEstimator=="V0"){
1145  centrality = centrV0M;
1146  }
1147  else if(sCentrEstimator=="CL1"){
1148  centrality = centrCL1;
1149  }
1150  else if(sCentrEstimator=="V0C"){
1151  centrality = fMultSelection->GetMultiplicityPercentile("V0C");
1152  }
1153  else if(sCentrEstimator=="V0A"){
1154  centrality = fMultSelection->GetMultiplicityPercentile("V0A");
1155  }
1156  else if(sCentrEstimator=="TRK"){
1157  centrality = fMultSelection->GetMultiplicityPercentile("TRK");
1158  } */
1159 
1160 
1161 
1162 
1164  return;
1165  }
1166 
1167  fHistEventCount->Fill(stepCount); //4
1168  stepCount++;
1169 
1170  fCentDistBefore->Fill(centrality);
1171 
1172 
1173 
1174 
1175 
1176 
1177  Int_t ntracks=fAOD->GetNumberOfTracks();
1178  if(ntracks<2) return; // Check this cut....!!!
1179 
1180  fHistEventCount->Fill(stepCount); //5
1181  stepCount++;
1182 
1183 
1184 
1185 
1186 
1187 
1188  Int_t cent10bin = -1;
1189  Int_t cIndex = -1;
1190  //cent10bin = GetCentralityScaled0to10(centrality); //Centrality in 0-10 scale
1191 
1192  if(centrality<5.0) {
1193  cent10bin = 0;
1194  }
1195  else if(centrality>=5.0 && centrality<10){
1196  cent10bin = 1;
1197  }
1198  else if(centrality>=10.0) {
1199  cent10bin = abs(centrality/10.0)+1;
1200  }
1201 
1202  cIndex = cent10bin;
1203 
1204 //Centrality array index for NUA correcion
1205  Int_t cForNUA = 0;
1206 
1207  if(centrality<5.0) {
1208  cForNUA = 0;
1209  }
1210  else if(centrality>=5.0 && centrality<10){
1211  cForNUA = 1; // 1=5-10,
1212  }
1213  else if(centrality>=10.0 && centrality<20) {
1214  cForNUA = 2; // 2 = 10-20,
1215  }
1216  else if(centrality>=20 && centrality<40){
1217  cForNUA = 3; // 3=20-40
1218  }
1219  else if(centrality>=40){
1220  cForNUA = 4; // 4=40-90
1221  }
1222 
1223 
1224 
1225 
1226 
1227 
1228  //---------- Magnetic field --------
1229  Double_t fMagField = fAOD->GetMagneticField();
1230 
1231  const Int_t QAindex = (fMagField > 0) ? 1 : 0;
1232  //---------------------------------
1233 
1234 
1235 
1236 
1237 
1238 
1239  //Load NUA and V0M correction map run by run:
1240  Int_t runNumber = fAOD->GetRunNumber();
1241 
1242  if(runNumber!=fOldRunNum) {
1243 
1244  GetNUACorrectionHist(runNumber);
1245 
1246  if(bV0MGainCorr) {
1247  GetV0MCorrectionHist(runNumber);
1248  }
1249 
1250  fOldRunNum = runNumber;
1251  }
1252  //------------------------------------------
1253 
1254 
1255 
1256 
1257 
1258 
1259 
1260  //----- Event Plane variables:-------
1261  Double_t PsiNV0A = 0;
1262  Double_t PsiNV0C = 0;
1263 
1264  Double_t PsiNTPCA = 0; // eta <0
1265  Double_t PsiNTPCC = 0; // eta >0
1266  Double_t PsiNTPCF = 0; // Full TPC
1267 
1268  Double_t sumTPCQn2x[3] = {0,0,0}; //[0]= eta<0; [1]= eta>0; [2]= -0.8 < eta < 0.8
1269  Double_t sumTPCQn2y[3] = {0,0,0};
1270  //Double_t sumTPCQn3x[3] = {0,0,0}; //[0]= eta<0; [1]= eta>0; [2]= -0.8 < eta < 0.8
1271  //Double_t sumTPCQn3y[3] = {0,0,0};
1272  //Double_t sumTPCQn4x[3] = {0,0,0}; //[0]= eta<0; [1]= eta>0; [2]= -0.8 < eta < 0.8
1273  //Double_t sumTPCQn4y[3] = {0,0,0};
1274  //------------------------------------
1275 
1276 
1277 
1278 
1279 
1280 
1281 
1282  //Variables for MC tracking correction
1283  Int_t ptBinMC = 1;
1284  Int_t iBinNUA = 1;
1285  Double_t ptWgtMC = 1.0;
1286  Double_t WgtNUA = 1.0;
1287  Double_t ptTrk = 0.1;
1288  Double_t dEdx = 0.0;
1289  Double_t Chi2Trk = 0.0;
1290 
1291 
1292  //-------------- Track loop for outlier and PileUp cut -------------------
1293 
1294  //---------------- a dobrin --------------
1295 
1296  Bool_t bIsPileup=kFALSE;
1297 
1298  Int_t isPileup = fAOD->IsPileupFromSPD(3);
1299 
1300  if(isPileup != 0) {
1301  fHistPileUpCount->Fill(0.5);
1302  bIsPileup=kTRUE;
1303  }
1304  else if(PileUpMultiVertex(fAOD)) {
1305  fHistPileUpCount->Fill(1.5);
1306  bIsPileup=kTRUE;
1307  }
1308  else if(((AliAODHeader*)fAOD->GetHeader())->GetRefMultiplicityComb08() < 0) {
1309  fHistPileUpCount->Fill(2.5);
1310  bIsPileup=kTRUE;
1311  }
1312  else if(fAOD->IsIncompleteDAQ()) {
1313  fHistPileUpCount->Fill(3.5);
1314  bIsPileup=kTRUE;
1315  }
1316  else if(fabs(centrV0M-centrCL1)> 5.0) {//default: 7.5
1317 //else if(fabs(centrV0M-centrCL1)> 7.5) {//default: 7.5
1318  fHistPileUpCount->Fill(4.5);
1319  bIsPileup=kTRUE;
1320  }
1321 
1322  // check vertex consistency
1323  const AliAODVertex* vtTrc = fAOD->GetPrimaryVertex();
1324  const AliAODVertex* vtSPD = fAOD->GetPrimaryVertexSPD();
1325 
1326  if(vtTrc->GetNContributors() < 2 || vtSPD->GetNContributors()<1) {
1327  fHistPileUpCount->Fill(5.5);
1328  bIsPileup=kTRUE;
1329  }
1330 
1331  double covTrc[6], covSPD[6];
1332  vtTrc->GetCovarianceMatrix(covTrc);
1333  vtSPD->GetCovarianceMatrix(covSPD);
1334 
1335  double dz = vtTrc->GetZ() - vtSPD->GetZ();
1336 
1337  double errTot = TMath::Sqrt(covTrc[5]+covSPD[5]);
1338  double errTrc = TMath::Sqrt(covTrc[5]);
1339  double nsigTot = dz/errTot;
1340  double nsigTrc = dz/errTrc;
1341 
1342  if(TMath::Abs(dz)>0.2 || TMath::Abs(nsigTot)>10 || TMath::Abs(nsigTrc)>20) {
1343  fHistPileUpCount->Fill(6.5);
1344  bIsPileup=kTRUE;
1345  }
1346 
1347 
1348  /*
1349  Float_t nSigPionTPC[40000] = {0.,};
1350  Float_t nSigKaonTPC[40000] = {0.,};
1351  Float_t nSigProtonTPC[40000] = {0.,};
1352  Float_t nSigPionTOF[40000] = {0.,};
1353  Float_t nSigKaonTOF[40000] = {0.,};
1354  Float_t nSigProtonTOF[40000] = {0.,};
1355  */
1356 
1357  if(ntracks > 40000) return; //Dont break segment for higher tracks:
1358 
1359 
1360  /*
1361  std::vector<float> nSigPionTPC;
1362  std::vector<float> nSigKaonTPC;
1363  std::vector<float> nSigProtonTPC;
1364  std::vector<float> nSigPionTOF;
1365  std::vector<float> nSigKaonTOF;
1366  std::vector<float> nSigProtonTOF;
1367 
1368  nSigPionTPC.reserve(10000);
1369  nSigKaonTPC.reserve(10000);
1370  nSigProtonTPC.reserve(10000);
1371  nSigPionTOF.reserve(10000);
1372  nSigKaonTOF.reserve(10000);
1373  nSigProtonTOF.reserve(10000); */
1374 
1375 
1376 
1377 
1378  Float_t multTPC = 0; // tpc mult estimate
1379  //Float_t RefMultRaw = 0; // tpc mult estimate
1380  //Float_t RefMultCorr = 0; // tpc mult estimate
1381  Float_t RefMultRawFB = 0;
1382  Float_t RefMultCorrFB= 0;
1383 
1384  Float_t multTPCAll = 0; // tpc mult estimate
1385  Float_t multGlobal = 0; // global multiplicity
1386 
1387  Int_t multEtaNeg, multEtaPos, multEtaFull;
1388  Double_t SumWEtaNeg, SumWEtaPos, SumWEtaFull;
1389 
1390  Int_t ChTrk;
1391  Int_t nClustTPC;
1392 
1393  Double_t etaTrk, phiTrk; // Never define eta as float, always double.!!
1394 
1395  //Double_t fMaxPtEP = 5.0;
1396  //Double_t fMinPtEP = 0.2;
1397  Double_t fMaxEtaEP = 0.8;
1398  Double_t fMinEtaEP = -0.8;
1399 
1400  multEtaNeg = 0;
1401  multEtaPos = 0;
1402  multEtaFull= 0;
1403  SumWEtaNeg = 0;
1404  SumWEtaPos = 0;
1405  SumWEtaFull= 0;
1406 
1407 
1408 
1409  for(Int_t iTrack = 0; iTrack < ntracks; iTrack++) { //-------------------------
1410 
1411  AliAODTrack* AODtrack =dynamic_cast<AliAODTrack*>(fVevent->GetTrack(iTrack));
1412  if(!AODtrack) continue;
1413 
1414 
1415  //------ Store PID info in array to use later:---------
1416  //Vector method: // taking same time as arrays
1417  /*
1418  nSigPionTPC.push_back(fPIDResponse->NumberOfSigmasTPC(AODtrack, AliPID::kPion));
1419  nSigKaonTPC.push_back(fPIDResponse->NumberOfSigmasTPC( AODtrack, AliPID::kKaon));
1420  nSigProtonTPC.push_back(fPIDResponse->NumberOfSigmasTPC(AODtrack, AliPID::kProton));
1421  nSigPionTOF.push_back(fPIDResponse->NumberOfSigmasTOF(AODtrack, AliPID::kPion));
1422  nSigKaonTOF.push_back(fPIDResponse->NumberOfSigmasTOF( AODtrack, AliPID::kKaon));
1423  nSigProtonTOF.push_back(fPIDResponse->NumberOfSigmasTOF(AODtrack, AliPID::kProton)); */
1424 
1425  //I used these Arrays for PID, now commented as not enough statistics: 28/09/18
1426  /*
1427  nSigPionTPC[iTrack] = fPIDResponse->NumberOfSigmasTPC(AODtrack, AliPID::kPion);
1428  nSigKaonTPC[iTrack] = fPIDResponse->NumberOfSigmasTPC(AODtrack, AliPID::kKaon);
1429  nSigProtonTPC[iTrack] = fPIDResponse->NumberOfSigmasTPC(AODtrack, AliPID::kProton);
1430  nSigPionTOF[iTrack] = fPIDResponse->NumberOfSigmasTOF(AODtrack, AliPID::kPion);
1431  nSigKaonTOF[iTrack] = fPIDResponse->NumberOfSigmasTOF(AODtrack, AliPID::kKaon);
1432  nSigProtonTOF[iTrack] = fPIDResponse->NumberOfSigmasTOF(AODtrack, AliPID::kProton);
1433  */
1434  //-----------------------------------------------------
1435 
1436 
1437  ptTrk = AODtrack->Pt();
1438  etaTrk = AODtrack->Eta();
1439  Chi2Trk = AODtrack->Chi2perNDF();
1440  nClustTPC = AODtrack->GetTPCNcls();
1441  //dEdx = AODtrack->GetDetPid()->GetTPCsignal(); // This one breaks the code if called before checking FilterBit.
1442 
1443 
1444  //Track cuts for POIs: Same to be used for EP calc.
1445  //if((dPt2 > fMaxPtCut) || (dPt2 < fMinPtCut) || (dEta2 > fMaxEtaCut) || (dEta2 < fMinEtaCut) || (dEdx2 < fdEdxMin) || (track2->GetTPCNcls() < fTPCclustMin) || (Chi2Trk2 < fTrkChi2Min) || (Chi2Trk2 > 4.0) || (track2->DCA() > fDCAxyMax)|| (track2->ZAtDCA() > fDCAzMax) || !(TMath::Abs(gCharge2)))
1446 
1447  //cuts for EP calculation:
1448  if(AODtrack->TestFilterBit(fFilterBit)){
1449 
1450  phiTrk = AODtrack->Phi();
1451  dEdx = AODtrack->GetDetPid()->GetTPCsignal();
1452  ChTrk = AODtrack->Charge();
1453 
1454  if(gPsiN > 0 && (ptTrk <= fMaxPtCut) && (ptTrk >= fMinPtCut) && (etaTrk <= fMaxEtaEP) && (etaTrk >= fMinEtaEP) && (dEdx >= fdEdxMin) && (nClustTPC >= fTPCclustMin) && (AODtrack->DCA() <= fDCAxyMax) && (AODtrack->ZAtDCA() <= fDCAzMax) && (Chi2Trk >= fTrkChi2Min) && (Chi2Trk <= 4.0) && TMath::Abs(ChTrk))
1455  {
1456 
1457  ptWgtMC = 1.0;
1458 
1459  if(fFB_Efficiency_Cent[cent10bin]){
1460  ptBinMC = fFB_Efficiency_Cent[cent10bin]->FindBin(ptTrk); //Charge independent MC correction atm.
1461  ptWgtMC = 1.0/fFB_Efficiency_Cent[cent10bin]->GetBinContent(ptBinMC);
1462  }
1463 
1464  RefMultRawFB++;
1465  RefMultCorrFB += ptWgtMC;
1466 
1467 
1468  //------> Get NUA weights for EP <----------
1469  WgtNUA = 1.0;
1470 
1471  if(ChTrk>0){
1472  if(fHCorrectNUApos[cForNUA]){
1473  iBinNUA = fHCorrectNUApos[cForNUA]->FindBin(pVtxZ,phiTrk,etaTrk);
1474  WgtNUA = fHCorrectNUApos[cForNUA]->GetBinContent(iBinNUA);
1475  }
1476  //else{ WgtNUA = 1.0; }
1477  }
1478  else{
1479  if(fHCorrectNUAneg[cForNUA]){
1480  iBinNUA = fHCorrectNUAneg[cForNUA]->FindBin(pVtxZ,phiTrk,etaTrk);
1481  WgtNUA = fHCorrectNUAneg[cForNUA]->GetBinContent(iBinNUA);
1482  }
1483  //else{ WgtNUA = 1.0; }
1484  }
1485 
1486  if(etaTrk < -0.05){
1487  sumTPCQn2x[0] += WgtNUA*TMath::Cos(gPsiN*phiTrk);
1488  sumTPCQn2y[0] += WgtNUA*TMath::Sin(gPsiN*phiTrk);
1489  //sumTPCQn3x[0] += WgtNUA*TMath::Cos(3*phiTrk);
1490  //sumTPCQn3y[0] += WgtNUA*TMath::Sin(3*phiTrk);
1491  //sumTPCQn4x[0] += WgtNUA*TMath::Cos(4*phiTrk);
1492  //sumTPCQn4y[0] += WgtNUA*TMath::Sin(4*phiTrk);
1493  multEtaNeg++;
1494  SumWEtaNeg += WgtNUA;
1495  }
1496  else if(etaTrk > 0.05){
1497  sumTPCQn2x[1] += WgtNUA*TMath::Cos(gPsiN*phiTrk);
1498  sumTPCQn2y[1] += WgtNUA*TMath::Sin(gPsiN*phiTrk);
1499  //sumTPCQn3x[1] += WgtNUA*TMath::Cos(3*phiTrk);
1500  //sumTPCQn3y[1] += WgtNUA*TMath::Sin(3*phiTrk);
1501  //sumTPCQn4x[1] += WgtNUA*TMath::Cos(4*phiTrk);
1502  //sumTPCQn4y[1] += WgtNUA*TMath::Sin(4*phiTrk);
1503  multEtaPos++;
1504  SumWEtaPos += WgtNUA;
1505  }
1506  sumTPCQn2x[2] += WgtNUA*TMath::Cos(gPsiN*phiTrk);
1507  sumTPCQn2y[2] += WgtNUA*TMath::Sin(gPsiN*phiTrk);
1508  //sumTPCQn3x[2] += WgtNUA*TMath::Cos(3*phiTrk);
1509  //sumTPCQn3y[2] += WgtNUA*TMath::Sin(3*phiTrk);
1510  //sumTPCQn4x[2] += WgtNUA*TMath::Cos(4*phiTrk);
1511  //sumTPCQn4y[2] += WgtNUA*TMath::Sin(4*phiTrk);
1512  multEtaFull++;
1513  SumWEtaFull += WgtNUA;
1514  }// track cuts
1515  } // AOD fiter bit
1516 
1517 
1518  if(AODtrack->TestFilterBit(128)) multTPCAll++; // A. Dobrin TPC vs ESD PileUp Cut.
1519 
1520  if(!(AODtrack->TestFilterBit(1))) continue; // cuts for Outlier as in FlowEvent Task
1521 
1522  dEdx = AODtrack->GetDetPid()->GetTPCsignal();
1523 
1524  if((ptTrk < 0.2) || (ptTrk > 5.0) || (TMath::Abs(etaTrk) > 0.8) || (nClustTPC < 70) || (dEdx < 10.0) || (Chi2Trk < 0.1)) continue;
1525 
1526  if(AODtrack->GetDetPid() && Chi2Trk > 0.2) multTPC++;
1527 
1528  if(!AODtrack->TestFilterBit(16) || AODtrack->Chi2perNDF() < 0.1) continue;
1529  Double_t b[2] = {-99., -99.};
1530  Double_t bCov[3] = {-99., -99., -99.};
1531  AliAODTrack copy(*AODtrack);
1532  if(copy.PropagateToDCA(fVevent->GetPrimaryVertex(), fVevent->GetMagneticField(), 100., b, bCov) && TMath::Abs(b[0]) < 0.3 && TMath::Abs(b[1]) < 0.3) multGlobal++;
1533 
1534  }//--- track loop outlier/PileUp ----
1535 
1536 
1537 
1538 
1539 
1540 
1541  Int_t multEsd = ((AliAODHeader*)fAOD->GetHeader())->GetNumberOfESDTracks();
1542  Float_t multESDTPCDiff = (Float_t) multEsd - fPileUpSlopeParm*multTPCAll;
1543 
1544  //cout<<" Info:UserExec() called ... I am after PU cut... event = "<<fEventCount<<" \n";
1545 
1546  if(multESDTPCDiff > fPileUpConstParm) {
1547  fHistPileUpCount->Fill(7.5);
1548  bIsPileup=kTRUE;
1549  }
1550  else if(bIsPileup==kFALSE) {
1551 
1552  if(!fMultSelection->GetThisEventIsNotPileup()){
1553  fHistMultSelPUCount->Fill(0.5);
1554  bIsPileup=kTRUE;
1555  }
1556  if(!fMultSelection->GetThisEventIsNotPileupMV()){
1557  fHistMultSelPUCount->Fill(1.5);
1558  bIsPileup=kTRUE;
1559  }
1560  if(!fMultSelection->GetThisEventIsNotPileupInMultBins()){
1561  fHistMultSelPUCount->Fill(2.5);
1562  bIsPileup=kTRUE;
1563  }
1564  if(!fMultSelection->GetThisEventHasNoInconsistentVertices()){
1565  fHistMultSelPUCount->Fill(2.5);
1566  bIsPileup=kTRUE;
1567  }
1568  if(!fMultSelection->GetThisEventPassesTrackletVsCluster()){
1569  fHistMultSelPUCount->Fill(2.5);
1570  bIsPileup=kTRUE;
1571  }
1572  if(!fMultSelection->GetThisEventIsNotIncompleteDAQ()){
1573  fHistMultSelPUCount->Fill(2.5);
1574  bIsPileup=kTRUE;
1575  }
1576  if(!fMultSelection->GetThisEventHasGoodVertex2016()){
1577  fHistMultSelPUCount->Fill(2.5);
1578  bIsPileup=kTRUE;
1579  }
1580  if(bIsPileup) fHistPileUpCount->Fill(9.5);
1581  }
1582  //-----------------------------------------------------------------
1583 
1584 
1585 
1586 
1587 
1588 
1589 
1590 
1591 
1592  fHistTPCVsESDTrkBefore->Fill(multTPCAll,multEsd); //A. Dobrin
1593  fHistTPCvsGlobalMultBefore->Fill(multGlobal,multTPC);
1594 
1595 
1596  Bool_t bIsOutLier=kFALSE;
1597 
1598  if(multTPC < (-20.0+1.15*multGlobal) || multTPC > (200.+1.45*multGlobal)) { bIsOutLier = kTRUE;}
1599 
1600  fHistEventCount->Fill(stepCount); //6
1601  stepCount++;
1602 
1603 
1604  fHistTPConlyVsCL1Before->Fill(centrCL1,multTPCAll);
1605  fHistTPConlyVsV0MBefore->Fill(centrV0M,multTPCAll);
1606  fHistGlobalVsV0MBefore->Fill(centrV0M, multGlobal);
1607 
1608 
1609  //if bSkipPileUpCut is kTRUE then don't apply PileUp removal.
1610  if(!bSkipPileUpCut && bIsOutLier) return; //outlier TPC vs Global
1611 
1612  fHistTPCvsGlobalMultAfter->Fill(multGlobal,multTPC);
1613 
1614  fHistEventCount->Fill(stepCount); //7
1615  stepCount++;
1616 
1617 
1618  if(!bSkipPileUpCut && bIsPileup) return; //PileUp A. Dobrin
1619 
1620  fHistTPCVsESDTrkAfter->Fill(multTPCAll,multEsd);
1621 
1622  fHistEventCount->Fill(stepCount); //8
1623  stepCount++;
1624 
1625  //cout<<"After PU cut multTPC = "<<multTPC<<" multGlobal = "<<multGlobal<<endl;
1626 
1627 
1628 
1629 
1630 
1631  /*
1632  Int_t icentBin = centrality;
1633  icentBin++;
1634 
1635  Float_t TPCmultLowLimit = hCentvsTPCmultCuts->GetBinContent(icentBin,1);
1636  Float_t TPCmultHighLimit = hCentvsTPCmultCuts->GetBinContent(icentBin,1);
1637 
1638  TPCmultLowLimit -= 5.0 * hCentvsTPCmultCuts->GetBinContent(icentBin,2); //mean - 5sigma
1639  TPCmultHighLimit += 5.0 * hCentvsTPCmultCuts->GetBinContent(icentBin,2); //mean + 5sigma
1640  //std::cout<<" Cent = "<<centrality<<"\t icent = "<<icentBin<<" low = "<<TPCmultLowLimit<<"\t high = "<<TPCmultHighLimit<<std::endl;
1641 
1642  //centrality outlier
1643  if(!bSkipPileUpCut){ if(multTPC<TPCmultLowLimit || multTPC>TPCmultHighLimit) return; }
1644  */
1645 
1646 
1647 
1648 
1649 
1650  fHistEventCount->Fill(stepCount); //9
1651  stepCount++;
1652 
1653 
1654 
1655  fHistTPConlyVsCL1After->Fill(centrCL1,multTPCAll);
1656  fHistTPConlyVsV0MAfter->Fill(centrV0M,multTPCAll);
1657  fHistGlobalVsV0MAfter->Fill(centrV0M, multGlobal);
1658 
1659 
1660  // MC corrected Refmult:
1661  fHistRawVsCorrMultFB->Fill(RefMultRawFB,RefMultCorrFB); // FB set by AddTask..
1662 
1663 
1664  Float_t EvtCent = centrality;
1665 
1666 
1667 
1668 
1669  if(gPsiN > 0){
1670  if(multEtaNeg<2 || multEtaPos<2) return; //Minimum 2 tracks in each eta
1671  }
1672 
1673  fHistEventCount->Fill(stepCount); //10
1674  stepCount++;
1675 
1676  //--------------------------------------------------------
1677 
1678 
1679 
1680 
1681 
1682 
1683 
1684 
1685 
1686 
1687 
1688 
1689 
1690 
1691 
1692 
1693 
1694 
1695  //--------------- cent CL1 <= 90 cut ------------------
1696  Int_t icentV0Qn = centrCL1; // cent CL1 used for V0 calibration.
1697  icentV0Qn += 1;
1698 
1699  if(icentV0Qn>90) return;
1700 
1701 
1702 
1703  fCentDistVsVz->Fill(pVtxZ,centrality);
1704 
1705 
1706 
1707  fHistEventCount->Fill(stepCount); //11
1708  stepCount++;
1709 
1710 
1711 
1712 
1713  //=============== Get the ZDC data ====================
1714  const Int_t fNumberOfZDCChannels = 5;
1715 
1716  if(gPsiN>0) { // Fill ZDC info for gPsiN > 0
1717 
1718  AliAODZDC *aodZDC = fAOD->GetZDCData();
1719 
1720  if(!aodZDC) {
1721  cout<<"\n ********* Error: could not find ZDC data ************ \n "<<endl;
1722  }
1723 
1724  const Double_t *gZNATowerRawAOD = aodZDC->GetZNATowerEnergy();
1725  const Double_t *gZNCTowerRawAOD = aodZDC->GetZNCTowerEnergy();
1726  const Double_t *gZPATowerRawAOD = aodZDC->GetZPATowerEnergy();
1727  const Double_t *gZPCTowerRawAOD = aodZDC->GetZPCTowerEnergy();
1728 
1729  Double_t gZNATowerRaw[5] = {0.,0.,0.,0.,0.};
1730  Double_t gZNCTowerRaw[5] = {0.,0.,0.,0.,0.};
1731  Double_t gZPATowerRaw[5] = {0.,0.,0.,0.,0.};
1732  Double_t gZPCTowerRaw[5] = {0.,0.,0.,0.,0.};
1733 
1734  for(Int_t iChannel = 0; iChannel < fNumberOfZDCChannels; iChannel++) {
1735  gZNATowerRaw[iChannel] = gZNATowerRawAOD[iChannel];
1736  gZNCTowerRaw[iChannel] = gZNCTowerRawAOD[iChannel];
1737  gZPATowerRaw[iChannel] = gZPATowerRawAOD[iChannel];
1738  gZPCTowerRaw[iChannel] = gZPCTowerRawAOD[iChannel];
1739 
1740  fZPASignalPerChVsCent->Fill(centrCL1,iChannel,gZPATowerRaw[iChannel]);
1741  fZPCSignalPerChVsCent->Fill(centrCL1,iChannel,gZPCTowerRaw[iChannel]);
1742  fZNASignalPerChVsCent->Fill(centrCL1,iChannel,gZNATowerRaw[iChannel]);
1743  fZNCSignalPerChVsCent->Fill(centrCL1,iChannel,gZNCTowerRaw[iChannel]);
1744 
1745  }// ZDC channel loop
1746 
1747  }
1748 
1749 
1750 
1751  //-------- V0M info ---------------
1752 
1753  //do v0m recentering
1754  Double_t QxanCor = 0, QyanCor = 0;
1755  Double_t QxcnCor = 0, QycnCor = 0;
1756 
1757  Double_t Qxan3 = 0., Qyan3 = 0.;
1758  Double_t Qxcn3 = 0., Qycn3 = 0.;
1759  Double_t Qxan2 = 0., Qyan2 = 0.;
1760  Double_t Qxcn2 = 0., Qycn2 = 0.;
1761 
1762  Double_t phiV0;
1763  Float_t fMultv0 = 0;
1764  Float_t sumMa = 0;
1765  Float_t sumMc = 0;
1766 
1767 
1768  if(gPsiN>0) { // Calculate EP only for gPsiN > 0
1769 
1770  const AliAODVZERO *fAODV0 = fAOD->GetVZEROData();
1771 
1772  for(int iV0 = 0; iV0 < 64; iV0++) { //0-31 is V0C, 32-63 VOA
1773 
1774  fMultv0 = fAODV0->GetMultiplicity(iV0);
1775 
1776  fV0MultChVsRun->Fill(iV0+0.5,centrCL1,fMultv0);
1777 
1778  if(fHCorrectV0M){
1779  fMultv0 = fMultv0 * fHCorrectV0M->GetBinContent(iV0+1); // Gain Correction
1780  //cout<<"info: run = "<<runNumber<<" cent = "<<centrCL1<<"\t channel = "<<iV0<<" gain = "<<fHCorrectV0M->GetBinContent(iV0+1)<<endl;
1781  }
1782 
1783 
1784  //fV0MultChVsRun->Fill(iV0+0.5,runindex,fMultv0);
1785 
1786  phiV0 = TMath::PiOver4()*(0.5 + iV0 % 8);
1787 
1788  if(iV0 < 32){
1789  Qxcn2 += TMath::Cos(2*phiV0) * fMultv0;
1790  Qycn2 += TMath::Sin(2*phiV0) * fMultv0;
1791  Qxcn3 += TMath::Cos(3*phiV0) * fMultv0;
1792  Qycn3 += TMath::Sin(3*phiV0) * fMultv0;
1793  sumMc += fMultv0;
1794  }
1795  else if(iV0 >= 32){
1796  Qxan2 += TMath::Cos(2*phiV0) * fMultv0;
1797  Qyan2 += TMath::Sin(2*phiV0) * fMultv0;
1798  Qxan3 += TMath::Cos(3*phiV0) * fMultv0;
1799  Qyan3 += TMath::Sin(3*phiV0) * fMultv0;
1800  sumMa += fMultv0;
1801  }
1802  }//----- channel loop ----------
1803 
1804 
1805  if(gPsiN==3){
1806  QxanCor = Qxan3/sumMa; //3rd order event plane
1807  QyanCor = Qyan3/sumMa;
1808  QxcnCor = Qxcn3/sumMc;
1809  QycnCor = Qycn3/sumMc;
1810 
1811  if(fHAvgerageQnV0C && fHAvgerageQnV0A && icentV0Qn < 91){
1812  QxanCor -= fHAvgerageQnV0A->GetBinContent(icentV0Qn,3); //x = Cos
1813  QxcnCor -= fHAvgerageQnV0C->GetBinContent(icentV0Qn,3); //x = Cos
1814  QyanCor -= fHAvgerageQnV0A->GetBinContent(icentV0Qn,4); //y = Sin
1815  QycnCor -= fHAvgerageQnV0C->GetBinContent(icentV0Qn,4); //y = Sin
1816 
1817  fV0AQ3xVsCentRun->Fill(centrCL1,QxanCor);
1818  fV0AQ3yVsCentRun->Fill(centrCL1,QyanCor);
1819  fV0CQ3xVsCentRun->Fill(centrCL1,QxcnCor);
1820  fV0CQ3yVsCentRun->Fill(centrCL1,QycnCor);
1821  }
1822  //printf("\n .... I am using my own V0 gain correction for Psi3...\n");
1823  }
1824  else{
1825  QxanCor = Qxan2/sumMa; //2nd order Event plane
1826  QyanCor = Qyan2/sumMa;
1827  QxcnCor = Qxcn2/sumMc;
1828  QycnCor = Qycn2/sumMc;
1829 
1830  if(fHAvgerageQnV0C && fHAvgerageQnV0A && icentV0Qn < 91){
1831  QxanCor -= fHAvgerageQnV0A->GetBinContent(icentV0Qn,1); //x = Cos
1832  QxcnCor -= fHAvgerageQnV0C->GetBinContent(icentV0Qn,1); //x = Cos
1833  QyanCor -= fHAvgerageQnV0A->GetBinContent(icentV0Qn,2); //y = Sin
1834  QycnCor -= fHAvgerageQnV0C->GetBinContent(icentV0Qn,2); //y = Sin
1835 
1836  fV0AQ2xVsCentRun->Fill(centrCL1,QxanCor);
1837  fV0AQ2yVsCentRun->Fill(centrCL1,QyanCor);
1838  fV0CQ2xVsCentRun->Fill(centrCL1,QxcnCor);
1839  fV0CQ2yVsCentRun->Fill(centrCL1,QycnCor);
1840  }
1841  //printf("\n .... I am using my own V0 gain correction for Psi2...\n");
1842  }
1843 
1844  //------ For V0-Qn Recenter and Event plane: Uncorrectd ----------
1846  fV0CQ2xVsCentRun->Fill(centrCL1,Qxcn2/sumMc);
1847  fV0CQ2yVsCentRun->Fill(centrCL1,Qycn2/sumMc);
1848  fV0AQ2xVsCentRun->Fill(centrCL1,Qxan2/sumMa);
1849  fV0AQ2yVsCentRun->Fill(centrCL1,Qyan2/sumMa);
1850 
1851  fV0CQ3xVsCentRun->Fill(centrCL1,Qxcn3/sumMc);
1852  fV0CQ3yVsCentRun->Fill(centrCL1,Qycn3/sumMc);
1853  fV0AQ3xVsCentRun->Fill(centrCL1,Qxan3/sumMa);
1854  fV0AQ3yVsCentRun->Fill(centrCL1,Qyan3/sumMa);
1855  }
1856 
1857  if(gPsiN>2){
1858  PsiNV0C = 1.0/gPsiN*( TMath::ATan2(QycnCor,QxcnCor) + TMath::Pi() );
1859  PsiNV0A = 1.0/gPsiN*( TMath::ATan2(QyanCor,QxanCor) + TMath::Pi() );
1860  }
1861  else{
1862  PsiNV0C = 1.0/gPsiN*TMath::ATan2(QycnCor,QxcnCor) ;
1863  if(PsiNV0C<0.) PsiNV0C += 2*TMath::Pi()/gPsiN;
1864 
1865  PsiNV0A = 1.0/gPsiN*TMath::ATan2(QyanCor,QxanCor) ;
1866  if(PsiNV0A<0.) PsiNV0A += 2*TMath::Pi()/gPsiN;
1867  }
1868 
1869 
1870  fHV0CEventPlaneVsCent->Fill(EvtCent,PsiNV0C);
1871  fHV0AEventPlaneVsCent->Fill(EvtCent,PsiNV0A);
1872 
1873  } //--------------------- if(gPsiN>0) --------------------------
1874 
1875 
1876 
1877 
1878 
1879 
1880 
1881  fEventCount++;
1882 
1883 
1884 
1885 
1886 
1887 
1888 
1889 
1890 
1891  //---- Copies of TPC-Q vectors to remove track -by- track auto-correlation -----
1892 
1893  Double_t sumQxTPCneg;
1894  Double_t sumQyTPCneg;
1895  Double_t sumQxTPCpos;
1896  Double_t sumQyTPCpos;
1897  Double_t sumQxTPCneg2;
1898  Double_t sumQyTPCneg2;
1899  Double_t sumQxTPCpos2;
1900  Double_t sumQyTPCpos2;
1901  Double_t SumWgtNeg, SumWgtPos;
1902  Double_t SumWgtNeg2, SumWgtPos2;
1903 
1904 
1905 
1906  //--------- Track variable for PID/Charge studies ----------------
1907  Double_t PDGmassPion = 0.13957;
1908  Double_t PDGmassKaon = 0.49368;
1909  Double_t PDGmassProton = 0.93827;
1910 
1911  PDGmassProton *= PDGmassProton;
1912  PDGmassPion *= PDGmassPion;
1913  PDGmassKaon *= PDGmassKaon;
1914 
1915  Double_t dEdx1, dEdx2, dPhi1, dPhi2, dPt1, dPt2;
1916  Double_t dEta1, dEta2, ptw1, ptw2, deltaPhi;
1917  Double_t mom, w1NUA, w2NUA, WgtEP;
1918 
1919  /*
1920  Double_t nSigTOFpion, nSigTPCpion;
1921  Double_t nSigTOFkaon, nSigTPCkaon;
1922  Double_t nSigTOFproton,nSigTPCproton;
1923  Double_t nSigTOFpion2, nSigTPCpion2;
1924  Double_t nSigTOFkaon2, nSigTPCkaon2;
1925  Double_t nSigTOFproton2,nSigTPCproton2;
1926  */
1927 
1928  //Tof variables
1929  //Double_t length, tofTime, probMis, mass, beta;
1930  //Double_t c = TMath::C()*1.E-9; //bright light m/ns
1931  //Int_t TOFmatch=0;
1932 
1933 
1934 
1935  Int_t ptBin,gCharge1,gCharge2;
1936 
1937  //Double_t dcaXY, dcaZ ;
1938 
1939 
1940  //----------- Set the desired Harmonic ------------
1941  Int_t n = gN;
1942  Int_t m = gM;
1943  Int_t p =n+m;
1944  //------------------------------------------------
1945 
1946 
1947  Int_t skipPairHBT = 0;
1948 
1949 
1950 
1951 
1952 
1953  Double_t Chi2Trk1=0, Chi2Trk2=0;
1954 
1955  //Double_t ptwPion1,ptwKaon1,ptwProton1;
1956  //Double_t ptwPion2,ptwKaon2,ptwProton2;
1957  //Double_t wNUAPion1,wNUAKaon1,wNUAProton1;
1958  //Double_t wNUAPion2,wNUAKaon2,wNUAProton2;
1959 
1960  //Double_t WgtEPPion = 1.0;
1961  //Double_t WgtEPKaon = 1.0;
1962  //Double_t WgtEPProton = 1.0;
1963 
1964 
1965  Bool_t isPion1 = kFALSE;
1966  Bool_t isKaon1 = kFALSE;
1967  Bool_t isProton1 = kFALSE;
1968  Bool_t isPion2 = kFALSE;
1969  Bool_t isKaon2 = kFALSE;
1970  Bool_t isProton2 = kFALSE;
1971 
1972  Int_t multPOI1st = 0;
1973  Int_t multPOI2nd = 0;
1974 
1975 
1976  //const Int_t maxTrack = 40000;
1977 
1978 
1979 
1980 
1981  //Calling fPIDResponse in nested loop is CPU expensive.
1982  //Store nSigma values in a array:
1983 
1984 
1985  /*
1986  for(int itrack = 0; itrack < ntracks; itrack++) {
1987 
1988  AliAODTrack *trackForPID=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(itrack));
1989  if(!trackForPID) continue;
1990 
1991  // Array method:
1992  if(trackForPID){
1993  nSigPionTPC[itrack] = fPIDResponse->NumberOfSigmasTPC(trackForPID, AliPID::kPion);
1994  nSigKaonTPC[itrack] = fPIDResponse->NumberOfSigmasTPC(trackForPID, AliPID::kKaon);
1995  nSigProtonTPC[itrack] = fPIDResponse->NumberOfSigmasTPC(trackForPID, AliPID::kProton);
1996  nSigPionTOF[itrack] = fPIDResponse->NumberOfSigmasTOF(trackForPID, AliPID::kPion);
1997  nSigKaonTOF[itrack] = fPIDResponse->NumberOfSigmasTOF(trackForPID, AliPID::kKaon);
1998  nSigProtonTOF[itrack] = fPIDResponse->NumberOfSigmasTOF(trackForPID, AliPID::kProton);
1999  }
2000  else{
2001  nSigPionTPC[itrack] = -99;
2002  nSigKaonTPC[itrack] = -99;
2003  nSigProtonTPC[itrack] = -99;
2004  nSigPionTOF[itrack] = -99;
2005  nSigKaonTOF[itrack] = -99;
2006  nSigProtonTOF[itrack] = -99;
2007  }
2008 
2009  //if(itrack%10==0)
2010  //cout<< "nSig pi = " <<nSigPionTPC[itrack] << "nSig K = " << nSigKaonTPC[itrack] << "nSig p = " << nSigProtonTOF[itrack] << endl;
2011 
2012  // Vector method:
2013  //if(trackForPID){
2014  //nSigPionTPC.push_back(fPIDResponse->NumberOfSigmasTPC(trackForPID, AliPID::kPion));
2015  //nSigKaonTPC.push_back(fPIDResponse->NumberOfSigmasTPC(trackForPID, AliPID::kKaon));
2016  //nSigProtonTPC.push_back(fPIDResponse->NumberOfSigmasTPC(trackForPID, AliPID::kProton));
2017 
2018  //nSigPionTOF.push_back(fPIDResponse->NumberOfSigmasTOF(trackForPID, AliPID::kPion));
2019  //nSigKaonTOF.push_back(fPIDResponse->NumberOfSigmasTOF(trackForPID, AliPID::kKaon));
2020  //nSigProtonTOF.push_back(fPIDResponse->NumberOfSigmasTOF(trackForPID, AliPID::kProton));
2021  //}
2022  //else{
2023  //nSigPionTPC.push_back(-100.);
2024  //nSigKaonTPC.push_back(-100.);
2025  //nSigProtonTPC.push_back(-100.);
2026 
2027  //nSigPionTOF.push_back(-100.);
2028  //nSigKaonTOF.push_back(-100.);
2029  //nSigProtonTOF.push_back(-100.);
2030  //}
2031  }*/
2032 
2033 
2034 
2035 
2036 
2037 
2038  //Correct with Centrality Wgts : in second pass
2039  Int_t iCentBinWgt = (Int_t) centrality;
2040  iCentBinWgt += 1;
2041 
2042  Double_t fWgtCent = 1.0;
2043 
2044 
2045  if(fHCentWeightForRun){
2046  fWgtCent = fHCentWeightForRun->GetBinContent(iCentBinWgt);
2047  }
2048 
2049  //Fill Centrality for run-by-run: in first pass over data
2050  fCentDistAfter->Fill(centrality,fWgtCent);
2051 
2052 
2053 
2054 
2055 
2056 
2057  //#########################################
2058  // Track loop starting //
2059  //#########################################
2060 
2061 
2062 
2063 
2064  for(Int_t itrack = 0; itrack < ntracks; itrack++) {
2065 
2066  AliAODTrack *track=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(itrack));
2067  if(!track) continue;
2068 
2069  if(!track->TestFilterBit(fFilterBit)) continue;
2070 
2071 
2072 
2073  mom = track->P();
2074  dPt1 = track->Pt();
2075  dPhi1 = track->Phi();
2076  dEta1 = track->Eta();
2077  gCharge1 = track->Charge();
2078  dEdx1 = track->GetDetPid()->GetTPCsignal();
2079  Chi2Trk1 = track->Chi2perNDF();
2080 
2081  //dcaXY = track->DCA();
2082  //dcaZ = track->ZAtDCA();
2083 
2084  /* Turn off QAs
2085  //-------------- Check TOF status ------------------
2086  AliPIDResponse::EDetPidStatus status;
2087  status = fPIDResponse->CheckPIDStatus(AliPIDResponse::kTOF,track);
2088  TOFmatch = 0;
2089  if(status==AliPIDResponse::kDetPidOk){
2090  TOFmatch++;
2091  }
2092  probMis = fPIDResponse->GetTOFMismatchProbability(track);
2093  fHistTOFMatchCount->Fill(TOFmatch,probMis);
2094 
2095  mass = -9.9;
2096  beta = -0.9;
2097 
2098  if(TOFmatch>0 && probMis < 0.01) {
2099  //This conditions are called when detector status is checked above :
2100  //if((track->IsOn(AliAODTrack::kTOFin)) && (track->IsOn(AliAODTrack::kTIME)) && (track->IsOn(AliAODTrack::kTOFout))) {
2101  //if((track->IsOn(AliAODTrack::kITSin)) && (track->IsOn(AliAODTrack::kTOFpid))) { //Naghmeh used it
2102  tofTime = track->GetTOFsignal(); // in pico seconds
2103  length = track->GetIntegratedLength();
2104  tofTime = tofTime*1.E-3; // ns
2105  if (tofTime <= 0) tofTime = 9999;
2106  length = length*0.01; // in meters
2107  tofTime = tofTime*c;
2108  beta = length/tofTime;
2109  mass = mom*mom*(1./(beta*beta) - 1);
2110  }//------------ TOF signal -------------------------
2111 
2112 
2113  //QA histograms before applying track cuts:
2114  fHistEtaPtBefore->Fill(dEta1,dPt1);
2115  fHistTPCdEdxvsPBefore->Fill(mom*gCharge1,dEdx1);
2116  fHistTOFBetavsPBefore->Fill(mom*gCharge1,beta);
2117  fHistTOFMassvsPtBefore->Fill(dPt1*gCharge1,mass); */
2118 
2119 
2120 
2121 
2122  //-------- Apply Default track cuts for analysis: ---------
2123  if((dPt1 > fMaxPtCut) || (dPt1 < fMinPtCut) || (dEta1 > fMaxEtaCut) || (dEta1 < fMinEtaCut) || (dEdx1 < fdEdxMin) || (track->GetTPCNcls() < fTPCclustMin) || (Chi2Trk1 < fTrkChi2Min) || (Chi2Trk1 > 4.0) || (track->DCA() > fDCAxyMax) || (track->ZAtDCA() > fDCAzMax) || !(TMath::Abs(gCharge1)))
2124  continue;
2125 
2126  multPOI1st++;
2127 
2128 
2129 
2130  //--------------------- PID signals 1st track-------------------------
2131  //Vector/Array both works same way:
2132  /*
2133  nSigTPCpion = nSigPionTPC[itrack];
2134  nSigTPCkaon = nSigKaonTPC[itrack];
2135  nSigTPCproton = nSigProtonTPC[itrack];
2136  nSigTOFpion = nSigPionTOF[itrack];
2137  nSigTOFkaon = nSigKaonTOF[itrack];
2138  nSigTOFproton = nSigProtonTOF[itrack];
2139  */
2140 
2141 
2142  //if(itrack%100==0)
2143  //cout<<"Trk "<<itrack<<" pt1 = "<<dPt1<<"\tnSigPion = "<<nSigTPCpion<<"\tnSigKaon = "<<nSigTPCkaon <<"\tnSigprot = "<<nSigTPCproton<<endl;
2144 
2145 
2146  isPion1 = kFALSE;
2147  isKaon1 = kFALSE;
2148  isProton1 = kFALSE;
2149 
2150  /*
2151  //------> Pion
2152  if(dPt1<=0.6 && TMath::Abs(nSigTPCpion)<=3.0){
2153  isPion1 = kTRUE;
2154  }
2155  else if(dPt1>0.6 && dPt1<=2.0 && TMath::Abs(nSigTPCpion)<=3.0 && TMath::Abs(nSigTOFpion)<=2.0 ){
2156  isPion1 = kTRUE;
2157  }
2158  //------> Kaon
2159  if(dPt1<=0.45 && TMath::Abs(nSigTPCkaon)<=3.0){
2160  isKaon1 = kTRUE;
2161  }
2162  else if(dPt1>0.45 && dPt1<=2.0 && TMath::Abs(nSigTPCkaon)<=3.0 && TMath::Abs(nSigTOFkaon)<=2.0){
2163  isKaon1 = kTRUE;
2164  }
2165  //------> Proton
2166  if(dPt1<=0.8 && TMath::Abs(nSigTPCproton)<=3.0){
2167  isProton1 = kTRUE;
2168  if(dPt1<0.4) isProton1 = kFALSE; //Proton below 0.4 GeV is garbage
2169  }
2170  else if(dPt1>0.8 && dPt1<=2.0 && TMath::Abs(nSigTPCproton)<=3.0 && TMath::Abs(nSigTOFproton)<=3.0){ //0.8 to 3.5 was used earlier.
2171  isProton1 = kTRUE;
2172  } */
2173 
2174 
2175 
2176 
2177 
2178 
2179 
2180 
2181 
2182  //================= MC wgt and NUA wgt for PID =================
2183  /*
2184  ptwPion1 = 1.0;
2185  ptwKaon1 = 1.0;
2186  ptwProton1 = 1.0;
2187  wNUAPion1 = 1.0;
2188  wNUAKaon1 = 1.0;
2189  wNUAProton1 = 1.0;
2190 
2191  //------ get MC weight and NUA for Pion track1 --------------
2192  if(isPion1){
2193  if(fFB_Efficiency_Pion_Cent[cent10bin]){ // <-------------------- !!!! WARNING: use Pion Efficiency file when available.
2194  ptBin = fFB_Efficiency_Pion_Cent[cent10bin]->FindBin(dPt1);
2195  ptwPion1 = 1.0/fFB_Efficiency_Pion_Cent[cent10bin]->GetBinContent(ptBin);
2196  }
2197  //else{ ptwPion1 = 1.0; }
2198 
2199  if(gCharge1>0){
2200  if(fHCorrectNUAposPion[cForNUA]){
2201  iBinNUA = fHCorrectNUAposPion[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2202  wNUAPion1 = fHCorrectNUAposPion[cForNUA]->GetBinContent(iBinNUA);
2203  }
2204  //else{ wNUAPion1 = 1.0; }
2205  }
2206  else{
2207  if(fHCorrectNUAnegPion[cForNUA]){
2208  iBinNUA = fHCorrectNUAnegPion[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2209  wNUAPion1 = fHCorrectNUAnegPion[cForNUA]->GetBinContent(iBinNUA);
2210  }
2211  //else{ wNUAPion1 = 1.0; }
2212  }
2213  }
2214 
2215  //------ get MC weight and NUA for Kaon track1 --------------
2216  if(isKaon1){
2217  if(fFB_Efficiency_Kaon_Cent[cent10bin]){ // <-------------------- !!!! WARNING: use Kaon Efficiency file when available.
2218  ptBin = fFB_Efficiency_Kaon_Cent[cent10bin]->FindBin(dPt1);
2219  ptwKaon1 = 1.0/fFB_Efficiency_Kaon_Cent[cent10bin]->GetBinContent(ptBin);
2220  }
2221  //else{ ptwKaon1 = 1.0; }
2222 
2223  if(gCharge1>0){
2224  if(fHCorrectNUAposKaon[cForNUA]){
2225  iBinNUA = fHCorrectNUAposKaon[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2226  wNUAKaon1 = fHCorrectNUAposKaon[cForNUA]->GetBinContent(iBinNUA);
2227  }
2228  //else{ wNUAKaon1 = 1.0; }
2229  }
2230  else{
2231  if(fHCorrectNUAnegKaon[cForNUA]){
2232  iBinNUA = fHCorrectNUAnegKaon[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2233  wNUAKaon1 = fHCorrectNUAnegKaon[cForNUA]->GetBinContent(iBinNUA);
2234  }
2235  //else{ wNUAKaon1 = 1.0; }
2236  }
2237  }
2238 
2239  //------ get MC weight and NUA for Proton track1 --------------
2240  if(isProton1){
2241  if(gCharge1>0){
2242  if(fFB_Efficiency_Proton_Pos_Cent[cent10bin]){ // <-------------------- !!!! WARNING: use Proton Efficiency file when available.
2243  ptBin = fFB_Efficiency_Proton_Pos_Cent[cent10bin]->FindBin(dPt1);
2244  ptwProton1 = 1.0/fFB_Efficiency_Proton_Pos_Cent[cent10bin]->GetBinContent(ptBin);
2245  }
2246  if(fHCorrectNUAposProton[cForNUA]){
2247  iBinNUA = fHCorrectNUAposProton[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2248  wNUAProton1 = fHCorrectNUAposProton[cForNUA]->GetBinContent(iBinNUA);
2249  }
2250  }
2251  else{
2252  if(fFB_Efficiency_Proton_Neg_Cent[cent10bin]){ // <-------------------- !!!! WARNING: use Proton Efficiency file when available.
2253  ptBin = fFB_Efficiency_Proton_Neg_Cent[cent10bin]->FindBin(dPt1);
2254  ptwProton1 = 1.0/fFB_Efficiency_Proton_Neg_Cent[cent10bin]->GetBinContent(ptBin);
2255  }
2256  if(fHCorrectNUAnegProton[cForNUA]){
2257  iBinNUA = fHCorrectNUAnegProton[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2258  wNUAProton1 = fHCorrectNUAnegProton[cForNUA]->GetBinContent(iBinNUA);
2259  }
2260  }
2261  //else{ ptwProton1 = 1.0; }
2262  }
2263  */
2264  //=========================== X ===============================
2265 
2266 
2267 
2268  //------ get MC weight and NUA for Charged track 1--------------
2269  ptw1 = 1.0;
2270  w1NUA = 1.0;
2271 
2272  if(fFB_Efficiency_Cent[cent10bin]){
2273  ptBin = fFB_Efficiency_Cent[cent10bin]->FindBin(dPt1);
2274  ptw1 = 1.0/fFB_Efficiency_Cent[cent10bin]->GetBinContent(ptBin);
2275  }
2276  //else{ ptw1 = 1.0; }
2277 
2278  if(gCharge1>0){
2279  if(fHCorrectNUApos[cForNUA]){
2280  iBinNUA = fHCorrectNUApos[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2281  w1NUA = fHCorrectNUApos[cForNUA]->GetBinContent(iBinNUA);
2282  }
2283  //else{ w1NUA = 1.0; }
2284  }
2285  else{
2286  if(fHCorrectNUAneg[cForNUA]){
2287  iBinNUA = fHCorrectNUAneg[cForNUA]->FindBin(pVtxZ,dPhi1,dEta1);
2288  w1NUA = fHCorrectNUAneg[cForNUA]->GetBinContent(iBinNUA);
2289  }
2290  //else{ w1NUA = 1.0; }
2291  }
2292 
2293  if(w1NUA==0) w1NUA = 1.0;
2294 
2295 
2296 
2297  //============= charged hadron analysis: ============
2298  sumQxTPCneg = sumTPCQn2x[0]; // first copy to remove 1st track fron Q vector (AutoCorrelation)
2299  sumQyTPCneg = sumTPCQn2y[0];
2300  sumQxTPCpos = sumTPCQn2x[1];
2301  sumQyTPCpos = sumTPCQn2y[1];
2302  SumWgtNeg = SumWEtaNeg;
2303  SumWgtPos = SumWEtaPos;
2304 
2305  //-------- Remove phi1 from EP calculation ----------
2306 
2307  if(dEta1 < -0.05){
2308  sumQxTPCneg -= w1NUA*TMath::Cos(gPsiN*dPhi1);
2309  sumQyTPCneg -= w1NUA*TMath::Sin(gPsiN*dPhi1); // [0] = eta <-0.05
2310  SumWgtNeg -= w1NUA;
2311  }
2312  else if(dEta1 > 0.05){
2313  sumQxTPCpos -= w1NUA*TMath::Cos(gPsiN*dPhi1);
2314  sumQyTPCpos -= w1NUA*TMath::Sin(gPsiN*dPhi1); // [1] = eta > 0.05
2315  SumWgtPos -= w1NUA;
2316  }
2317  //-----------------------------------------------------
2318 
2319 
2320 
2321 
2322 
2323 
2324 
2325 
2326 
2327 
2328 
2329 
2330  multPOI2nd = 0;
2331 
2332  //---2nd track loop (nested)---
2333  for(Int_t jtrack = 0; jtrack < ntracks; jtrack++) {
2334 
2335  if(jtrack==itrack) continue;
2336 
2337  if(n==0) continue;
2338 
2339  AliAODTrack *track2=dynamic_cast<AliAODTrack*>(fVevent->GetTrack(jtrack));
2340  if(!track2) continue;
2341 
2342  if(!(track2->TestFilterBit(fFilterBit))) continue;
2343 
2344  dPt2 = track2->Pt();
2345  dPhi2 = track2->Phi();
2346  dEta2 = track2->Eta();
2347  gCharge2= track2->Charge();
2348  dEdx2 = track2->GetDetPid()->GetTPCsignal();
2349  Chi2Trk2 = track2->Chi2perNDF();
2350 
2351  if((dPt2 > fMaxPtCut) || (dPt2 < fMinPtCut) || (dEta2 > fMaxEtaCut) || (dEta2 < fMinEtaCut) || (dEdx2 < fdEdxMin) || (track2->GetTPCNcls() < fTPCclustMin) || (Chi2Trk2 < fTrkChi2Min) || (Chi2Trk2 > 4.0) || (track2->DCA() > fDCAxyMax) || (track2->ZAtDCA() > fDCAzMax) || !(TMath::Abs(gCharge2)))
2352  continue;
2353 
2354  multPOI2nd++;
2355 
2356 
2357 
2358  /*
2359  if(track2->GetDetPid()->GetTPCsignal() < 10.0) continue;
2360  if(track2->GetTPCNcls() < 70) continue;
2361  if(track2->Chi2perNDF() < 0.2) continue;
2362  if(dPt2 < fMinPtCut || dPt2 > fMaxPtCut) continue;
2363  if(dEta2< fMinEtaCut || dEta2 > fMaxEtaCut) continue;
2364  if(!TMath::Abs(gCharge2)) continue;
2365  //if(!(track2->TestFilterBit(fFilterBit))) continue; */
2366  //-----------------------------------------------------
2367 
2368  //calling fPIDResponse is too costly for CPU
2369 
2370 
2371  isPion2 = kFALSE;
2372  isKaon2 = kFALSE;
2373  isProton2 = kFALSE;
2374 
2375  //--------------------- PID signals 2nd track-------------------------
2376  //Vector/Array both works same way:
2377  /*
2378  nSigTPCpion2 = nSigPionTPC[jtrack];
2379  nSigTPCkaon2 = nSigKaonTPC[jtrack];
2380  nSigTPCproton2 = nSigProtonTPC[jtrack];
2381  nSigTOFpion2 = nSigPionTOF[jtrack];
2382  nSigTOFkaon2 = nSigKaonTOF[jtrack];
2383  nSigTOFproton2 = nSigProtonTOF[jtrack];
2384 
2385  //------> Pion
2386  if(dPt2<=0.6 && TMath::Abs(nSigTPCpion2)<=3.0){
2387  isPion2 = kTRUE;
2388  }
2389  else if(dPt2>0.6 && dPt2<=2.0 && TMath::Abs(nSigTPCpion2)<=3.0 && TMath::Abs(nSigTOFpion2)<=2.0 ){
2390  isPion2 = kTRUE;
2391  }
2392  //------> Kaon
2393  if(dPt2<=0.45 && TMath::Abs(nSigTPCkaon2)<=3.0){
2394  isKaon2 = kTRUE;
2395  }
2396  else if(dPt2>0.45 && dPt2<=2.0 && TMath::Abs(nSigTPCkaon2)<=3.0 && TMath::Abs(nSigTOFkaon2)<=2.0){
2397  isKaon2 = kTRUE;
2398  }
2399  //------> Proton
2400  if(dPt2<=0.8 && TMath::Abs(nSigTPCproton2)<=3.0){
2401  isProton2 = kTRUE;
2402  if(dPt2<0.4) isProton2 = kFALSE; //Proton below 0.4 GeV is garbage
2403  }
2404  else if(dPt2>0.8 && dPt2<=2.0 && TMath::Abs(nSigTPCproton2)<=3.0 && TMath::Abs(nSigTOFproton2)<=3.0){ //0.8 to 3.5 was used earlier.
2405  isProton2 = kTRUE;
2406  }
2407  */
2408 
2409  //----------------------------------------------------------------
2410 
2411 
2412 
2413 
2414 
2415  //================= MC wgt and NUA wgt for PID =================
2416  /*
2417  WgtEPPion = 1.0;
2418  WgtEPKaon = 1.0;
2419  WgtEPProton = 1.0;
2420 
2421  ptwPion2 = 1.0;
2422  ptwKaon2 = 1.0;
2423  ptwProton2 = 1.0;
2424  wNUAPion2 = 1.0;
2425  wNUAKaon2 = 1.0;
2426  wNUAProton2 = 1.0;
2427 
2428  //------ get MC weight and NUA for Pion track2 --------------
2429  if(isPion2){
2430  if(fFB_Efficiency_Pion_Cent[cent10bin]){
2431  ptBin = fFB_Efficiency_Pion_Cent[cent10bin]->FindBin(dPt2);
2432  ptwPion2 = 1.0/fFB_Efficiency_Pion_Cent[cent10bin]->GetBinContent(ptBin);
2433  }
2434  //else{ ptwPion2 = 1.0; }
2435 
2436  if(gCharge2>0){
2437  if(fHCorrectNUAposPion[cForNUA]){
2438  iBinNUA = fHCorrectNUAposPion[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2439  wNUAPion2 = fHCorrectNUAposPion[cForNUA]->GetBinContent(iBinNUA);
2440  }
2441  //else{ wNUAPion2 = 1.0; }
2442  }
2443  else{
2444  if(fHCorrectNUAnegPion[cForNUA]){
2445  iBinNUA = fHCorrectNUAnegPion[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2446  wNUAPion2 = fHCorrectNUAnegPion[cForNUA]->GetBinContent(iBinNUA);
2447  }
2448  //else{ wNUAPion2 = 1.0; }
2449  }
2450 
2451  WgtEPPion = ptwPion1*ptwPion2*wNUAPion1*wNUAPion2;
2452  }
2453 
2454  //------ get MC weight and NUA for Kaon track2 --------------
2455  if(isKaon2){
2456  if(fFB_Efficiency_Kaon_Cent[cent10bin]){
2457  ptBin = fFB_Efficiency_Kaon_Cent[cent10bin]->FindBin(dPt2);
2458  ptwKaon2 = 1.0/fFB_Efficiency_Kaon_Cent[cent10bin]->GetBinContent(ptBin);
2459  }
2460  //else{ ptwKaon2 = 1.0; }
2461 
2462  if(gCharge2>0){
2463  if(fHCorrectNUAposKaon[cForNUA]){
2464  iBinNUA = fHCorrectNUAposKaon[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2465  wNUAKaon2 = fHCorrectNUAposKaon[cForNUA]->GetBinContent(iBinNUA);
2466  }
2467  //else{ wNUAKaon2 = 1.0; }
2468  }
2469  else{
2470  if(fHCorrectNUAnegKaon[cForNUA]){
2471  iBinNUA = fHCorrectNUAnegKaon[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2472  wNUAKaon2 = fHCorrectNUAnegKaon[cForNUA]->GetBinContent(iBinNUA);
2473  }
2474  //else{ wNUAKaon2 = 1.0; }
2475  }
2476 
2477  WgtEPKaon = ptwKaon1*ptwKaon2*wNUAKaon1*wNUAKaon2;
2478  }
2479 
2480  //------ get MC weight and NUA for Proton track2 --------------
2481  if(isProton2){
2482  if(gCharge2>0){
2483  if(fFB_Efficiency_Proton_Pos_Cent[cent10bin]){
2484  ptBin = fFB_Efficiency_Proton_Pos_Cent[cent10bin]->FindBin(dPt2);
2485  ptwProton2 = 1.0/fFB_Efficiency_Proton_Pos_Cent[cent10bin]->GetBinContent(ptBin);
2486  }
2487  if(fHCorrectNUAposProton[cForNUA]){
2488  iBinNUA = fHCorrectNUAposProton[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2489  wNUAProton2 = fHCorrectNUAposProton[cForNUA]->GetBinContent(iBinNUA);
2490  }
2491  }
2492  else{
2493  if(fFB_Efficiency_Proton_Neg_Cent[cent10bin]){
2494  ptBin = fFB_Efficiency_Proton_Neg_Cent[cent10bin]->FindBin(dPt2);
2495  ptwProton2 = 1.0/fFB_Efficiency_Proton_Neg_Cent[cent10bin]->GetBinContent(ptBin);
2496  }
2497  if(fHCorrectNUAnegProton[cForNUA]){
2498  iBinNUA = fHCorrectNUAnegProton[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2499  wNUAProton2 = fHCorrectNUAnegProton[cForNUA]->GetBinContent(iBinNUA);
2500  }
2501  }
2502  //else{ ptwProton2 = 1.0; }
2503  WgtEPProton = ptwProton1*ptwProton2*wNUAProton1*wNUAProton2;
2504  }
2505  */
2506  //========================== X ================================
2507 
2508 
2509 
2510 
2511 
2512  //------ get MC weight and NUA (Charged) for track 2--------------
2513  WgtEP = 1.0;
2514  ptw2 = 1.0;
2515  w2NUA = 1.0;
2516 
2517  if(fFB_Efficiency_Cent[cent10bin]){
2518  ptBin = fFB_Efficiency_Cent[cent10bin]->FindBin(dPt2);
2519  ptw2 = 1.0/fFB_Efficiency_Cent[cent10bin]->GetBinContent(ptBin);
2520  }
2521  //else{ ptw2 = 1.0; }
2522 
2523  if(gCharge2>0){
2524  if(fHCorrectNUApos[cForNUA]){
2525  iBinNUA = fHCorrectNUApos[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2526  w2NUA = fHCorrectNUApos[cForNUA]->GetBinContent(iBinNUA);
2527  }
2528  //else{ w2NUA = 1.0; }
2529  }
2530  else{
2531  if(fHCorrectNUAneg[cForNUA]){
2532  iBinNUA = fHCorrectNUAneg[cForNUA]->FindBin(pVtxZ,dPhi2,dEta2);
2533  w2NUA = fHCorrectNUAneg[cForNUA]->GetBinContent(iBinNUA);
2534  }
2535  //else{ w2NUA = 1.0; }
2536  }
2537 
2538  if(w2NUA==0) w2NUA = 1.0;
2539 
2540 
2541  //---------- Remove track2 from EP calculation ---------
2542  sumQxTPCneg2 = sumQxTPCneg; //second copy to remove 2nd track from EP
2543  sumQyTPCneg2 = sumQyTPCneg;
2544  sumQxTPCpos2 = sumQxTPCpos;
2545  sumQyTPCpos2 = sumQyTPCpos;
2546  SumWgtNeg2 = SumWgtNeg;
2547  SumWgtPos2 = SumWgtPos;
2548  //------------------------------------------------------
2549 
2550 
2551 
2552  if(dEta2 < -0.05){
2553  sumQxTPCneg2 -= w2NUA*TMath::Cos(gPsiN*dPhi2);
2554  sumQyTPCneg2 -= w2NUA*TMath::Sin(gPsiN*dPhi2); // [0] = eta <-0.05
2555  SumWgtNeg2 -= w2NUA;
2556  }
2557  else if(dEta2 > 0.05){
2558  sumQxTPCpos2 -= w2NUA*TMath::Cos(gPsiN*dPhi2);
2559  sumQyTPCpos2 -= w2NUA*TMath::Sin(gPsiN*dPhi2); // [1] = eta > 0.05
2560  SumWgtPos2 -= w2NUA;
2561  }
2562 
2563  if(SumWgtNeg2>0 && SumWgtPos2>0){
2564  sumQyTPCneg2 = sumQyTPCneg2/SumWgtNeg2;
2565  sumQxTPCneg2 = sumQxTPCneg2/SumWgtNeg2;
2566 
2567  sumQyTPCpos2 = sumQyTPCpos2/SumWgtPos2;
2568  sumQxTPCpos2 = sumQxTPCpos2/SumWgtPos2;
2569  }
2570 
2571  // track by track EP:
2572  if(gPsiN>0){
2573  PsiNTPCC = (1.0/gPsiN)*( TMath::ATan2(sumQyTPCneg2,sumQxTPCneg2) ) ; // negetive eta
2574  if(PsiNTPCC<0.) PsiNTPCC += 2*TMath::Pi()/gPsiN;
2575 
2576  PsiNTPCA = (1.0/gPsiN)*( TMath::ATan2(sumQyTPCpos2,sumQxTPCpos2) ) ; // positive eta
2577  if(PsiNTPCA<0.) PsiNTPCA += 2*TMath::Pi()/gPsiN;
2578  }
2579 
2580  //fHTPCAEventPlaneVsCent->Fill(EvtCent,PsiNTPCA);
2581  //fHTPCCEventPlaneVsCent->Fill(EvtCent,PsiNTPCC);
2582  //-----------------------------------------------------------
2583 
2584 
2585 
2586 
2587  // combined weight for EP:
2588  WgtEP = ptw1*ptw2*w1NUA*w2NUA;
2589 
2590 
2591 
2592  deltaPhi = dPhi1 - dPhi2; //for 2p correlator
2593 
2595 
2596  if(gCharge1!=gCharge2) {
2597  fHist_Corr3p_EP_Norm_PN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEP*fWgtCent);
2598  fHist_Corr3p_EP_Norm_PN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEP*fWgtCent);
2599  fHist_Corr3p_EP_Norm_PN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2600  fHist_Corr3p_EP_Norm_PN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2601  /*
2602  fHist_Corr3p_EP_Refm_PN[QAindex][0]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEP*fWgtCent);
2603  fHist_Corr3p_EP_Refm_PN[QAindex][1]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEP*fWgtCent);
2604  fHist_Corr3p_EP_Refm_PN[QAindex][2]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2605  fHist_Corr3p_EP_Refm_PN[QAindex][3]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2606  */
2607  //Differential
2608  if(cIndex<6){
2609  fHist_Corr3p_pTSum_EP_V0A_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2610  fHist_Corr3p_pTSum_EP_V0C_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2611  fHist_Corr3p_pTDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2612  fHist_Corr3p_pTDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2613  fHist_Corr3p_EtaDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2614  fHist_Corr3p_EtaDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2615  }
2616  //-------------> PID CME ---------------
2617  //Pion:
2618  if(isPion1 && isPion2){
2619  /*
2620  fHist_Corr3p_Pion_EP_Norm_PN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPPion*fWgtCent);
2621  fHist_Corr3p_Pion_EP_Norm_PN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPPion*fWgtCent);
2622  fHist_Corr3p_Pion_EP_Norm_PN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPPion*fWgtCent);
2623  fHist_Corr3p_Pion_EP_Norm_PN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPPion*fWgtCent);
2624 
2625  //Differential:
2626  if(cIndex<6){
2627  fHist_Corr3p_Pion_pTSum_EP_V0A_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2628  fHist_Corr3p_Pion_pTSum_EP_V0C_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2629  fHist_Corr3p_Pion_pTDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2630  fHist_Corr3p_Pion_pTDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2631  fHist_Corr3p_Pion_EtaDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2632  fHist_Corr3p_Pion_EtaDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2633  }
2634  */
2635  }
2636  //Kaon:
2637  if(isKaon1 && isKaon2){
2638  /*
2639  fHist_Corr3p_Kaon_EP_Norm_PN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPKaon*fWgtCent);
2640  fHist_Corr3p_Kaon_EP_Norm_PN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPKaon*fWgtCent);
2641  fHist_Corr3p_Kaon_EP_Norm_PN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPKaon*fWgtCent);
2642  fHist_Corr3p_Kaon_EP_Norm_PN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPKaon*fWgtCent);
2643 
2644  //Differential:
2645  if(cIndex<6){
2646  fHist_Corr3p_Kaon_pTSum_EP_V0A_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2647  fHist_Corr3p_Kaon_pTSum_EP_V0C_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2648  fHist_Corr3p_Kaon_pTDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2649  fHist_Corr3p_Kaon_pTDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2650  fHist_Corr3p_Kaon_EtaDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2651  fHist_Corr3p_Kaon_EtaDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2652  }
2653  */
2654  }
2655  //Proton:
2656  if(isProton1 && isProton2){
2657  //cout<<"#pair pt1 = "<<dPt1<<"\tpt2 = "<<dPt2<<"\tnSigp1 = "<<nSigTPCproton<<"\tnSigp2 = "<<nSigTPCproton2<<endl;
2658  /*
2659  fHist_Corr3p_Proton_EP_Norm_PN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPProton*fWgtCent);
2660  fHist_Corr3p_Proton_EP_Norm_PN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPProton*fWgtCent);
2661  fHist_Corr3p_Proton_EP_Norm_PN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPProton*fWgtCent);
2662  fHist_Corr3p_Proton_EP_Norm_PN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPProton*fWgtCent);
2663 
2664  //Differential:
2665  if(cIndex<6){
2666  fHist_Corr3p_Proton_pTSum_EP_V0A_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2667  fHist_Corr3p_Proton_pTSum_EP_V0C_PN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2668  fHist_Corr3p_Proton_pTDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2669  fHist_Corr3p_Proton_pTDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2670  fHist_Corr3p_Proton_EtaDiff_EP_V0A_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2671  fHist_Corr3p_Proton_EtaDiff_EP_V0C_PN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2672  }
2673  */
2674  }
2675  //------------------------------------
2676  }
2677 
2678 
2679 
2681 
2682 
2683 
2684  else if(gCharge1>0 && gCharge2>0 && skipPairHBT==0) {
2685  fHist_Corr3p_EP_Norm_PP[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEP*fWgtCent);
2686  fHist_Corr3p_EP_Norm_PP[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEP*fWgtCent);
2687  fHist_Corr3p_EP_Norm_PP[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2688  fHist_Corr3p_EP_Norm_PP[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2689  /*
2690  fHist_Corr3p_EP_Refm_PP[QAindex][0]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEP*fWgtCent);
2691  fHist_Corr3p_EP_Refm_PP[QAindex][1]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEP*fWgtCent);
2692  fHist_Corr3p_EP_Refm_PP[QAindex][2]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2693  fHist_Corr3p_EP_Refm_PP[QAindex][3]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2694  */
2695 
2696  //Differential:
2697  if(cIndex<6){
2698  fHist_Corr3p_pTSum_EP_V0A_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2699  fHist_Corr3p_pTSum_EP_V0C_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2700  fHist_Corr3p_pTDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2701  fHist_Corr3p_pTDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2702  fHist_Corr3p_EtaDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2703  fHist_Corr3p_EtaDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2704  }
2705  //-------------> PID CME ---------------
2706  //Pion:
2707  if(isPion1 && isPion2){
2708  /*
2709  fHist_Corr3p_Pion_EP_Norm_PP[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPPion*fWgtCent);
2710  fHist_Corr3p_Pion_EP_Norm_PP[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPPion*fWgtCent);
2711  fHist_Corr3p_Pion_EP_Norm_PP[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPPion*fWgtCent);
2712  fHist_Corr3p_Pion_EP_Norm_PP[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPPion*fWgtCent);
2713 
2714  //Differential:
2715  if(cIndex<6){
2716  fHist_Corr3p_Pion_pTSum_EP_V0A_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2717  fHist_Corr3p_Pion_pTSum_EP_V0C_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2718  fHist_Corr3p_Pion_pTDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2719  fHist_Corr3p_Pion_pTDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2720  fHist_Corr3p_Pion_EtaDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2721  fHist_Corr3p_Pion_EtaDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2722  }
2723  */
2724  }
2725  //Kaon:
2726  if(isKaon1 && isKaon2){
2727  /*
2728  fHist_Corr3p_Kaon_EP_Norm_PP[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPKaon*fWgtCent);
2729  fHist_Corr3p_Kaon_EP_Norm_PP[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPKaon*fWgtCent);
2730  fHist_Corr3p_Kaon_EP_Norm_PP[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPKaon*fWgtCent);
2731  fHist_Corr3p_Kaon_EP_Norm_PP[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPKaon*fWgtCent);
2732 
2733  //Differential:
2734  if(cIndex<6){
2735  fHist_Corr3p_Kaon_pTSum_EP_V0A_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2736  fHist_Corr3p_Kaon_pTSum_EP_V0C_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2737  fHist_Corr3p_Kaon_pTDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2738  fHist_Corr3p_Kaon_pTDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2739  fHist_Corr3p_Kaon_EtaDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2740  fHist_Corr3p_Kaon_EtaDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2741  }
2742  */
2743  }
2744  //Proton:
2745  if(isProton1 && isProton2){
2746  /*
2747  fHist_Corr3p_Proton_EP_Norm_PP[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPProton*fWgtCent);
2748  fHist_Corr3p_Proton_EP_Norm_PP[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPProton*fWgtCent);
2749  fHist_Corr3p_Proton_EP_Norm_PP[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPProton*fWgtCent);
2750  fHist_Corr3p_Proton_EP_Norm_PP[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPProton*fWgtCent);
2751 
2752  //Differential:
2753  if(cIndex<6){
2754  fHist_Corr3p_Proton_pTSum_EP_V0A_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2755  fHist_Corr3p_Proton_pTSum_EP_V0C_PP[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2756  fHist_Corr3p_Proton_pTDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2757  fHist_Corr3p_Proton_pTDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2758  fHist_Corr3p_Proton_EtaDiff_EP_V0A_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2759  fHist_Corr3p_Proton_EtaDiff_EP_V0C_PP[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2760  }
2761  */
2762  }
2763  //------------------------------------
2764  }
2765 
2766 
2767 
2768 
2769 
2771 
2772 
2773  else if(gCharge1<0 && gCharge2<0 && skipPairHBT==0){
2774  fHist_Corr3p_EP_Norm_NN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEP*fWgtCent);
2775  fHist_Corr3p_EP_Norm_NN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEP*fWgtCent);
2776  fHist_Corr3p_EP_Norm_NN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2777  fHist_Corr3p_EP_Norm_NN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2778  /*
2779  fHist_Corr3p_EP_Refm_NN[QAindex][0]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEP*fWgtCent);
2780  fHist_Corr3p_EP_Refm_NN[QAindex][1]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEP*fWgtCent);
2781  fHist_Corr3p_EP_Refm_NN[QAindex][2]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2782  fHist_Corr3p_EP_Refm_NN[QAindex][3]->Fill(RefMultCorrFB, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2783  */
2784  if(cIndex<6){
2785  fHist_Corr3p_pTSum_EP_V0A_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2786  fHist_Corr3p_pTSum_EP_V0C_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2787  fHist_Corr3p_pTDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2788  fHist_Corr3p_pTDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2789  fHist_Corr3p_EtaDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEP*fWgtCent);
2790  fHist_Corr3p_EtaDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEP*fWgtCent);
2791  }
2792  //-------------> PID CME ---------------
2793  //Pion:
2794  if(isPion1 && isPion2){
2795  /*
2796  fHist_Corr3p_Pion_EP_Norm_NN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPPion*fWgtCent);
2797  fHist_Corr3p_Pion_EP_Norm_NN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPPion*fWgtCent);
2798  fHist_Corr3p_Pion_EP_Norm_NN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPPion*fWgtCent);
2799  fHist_Corr3p_Pion_EP_Norm_NN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPPion*fWgtCent);
2800 
2801  //Differential:
2802  if(cIndex<6){
2803  fHist_Corr3p_Pion_pTSum_EP_V0A_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2804  fHist_Corr3p_Pion_pTSum_EP_V0C_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2805  fHist_Corr3p_Pion_pTDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2806  fHist_Corr3p_Pion_pTDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2807  fHist_Corr3p_Pion_EtaDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPPion*fWgtCent);
2808  fHist_Corr3p_Pion_EtaDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPPion*fWgtCent);
2809  }
2810  */
2811  }
2812  //Kaon:
2813  if(isKaon1 && isKaon2){
2814  /*
2815  fHist_Corr3p_Kaon_EP_Norm_NN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPKaon*fWgtCent);
2816  fHist_Corr3p_Kaon_EP_Norm_NN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPKaon*fWgtCent);
2817  fHist_Corr3p_Kaon_EP_Norm_NN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPKaon*fWgtCent);
2818  fHist_Corr3p_Kaon_EP_Norm_NN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPKaon*fWgtCent);
2819 
2820  //Differential:
2821  if(cIndex<6){
2822  fHist_Corr3p_Kaon_pTSum_EP_V0A_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2823  fHist_Corr3p_Kaon_pTSum_EP_V0C_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2824  fHist_Corr3p_Kaon_pTDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2825  fHist_Corr3p_Kaon_pTDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2826  fHist_Corr3p_Kaon_EtaDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPKaon*fWgtCent);
2827  fHist_Corr3p_Kaon_EtaDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPKaon*fWgtCent);
2828  }
2829  */
2830  }
2831  //Proton:
2832  if(isProton1 && isProton2){
2833  /*
2834  fHist_Corr3p_Proton_EP_Norm_NN[QAindex][0]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A), WgtEPProton*fWgtCent);
2835  fHist_Corr3p_Proton_EP_Norm_NN[QAindex][1]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C), WgtEPProton*fWgtCent);
2836  fHist_Corr3p_Proton_EP_Norm_NN[QAindex][2]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCA),WgtEPProton*fWgtCent);
2837  fHist_Corr3p_Proton_EP_Norm_NN[QAindex][3]->Fill(EvtCent, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNTPCC),WgtEPProton*fWgtCent);
2838 
2839  //Differential:
2840  if(cIndex<6){
2841  fHist_Corr3p_Proton_pTSum_EP_V0A_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2842  fHist_Corr3p_Proton_pTSum_EP_V0C_NN[QAindex][cIndex]->Fill((dPt1+dPt2)*0.5, TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2843  fHist_Corr3p_Proton_pTDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2844  fHist_Corr3p_Proton_pTDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dPt1-dPt2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2845  fHist_Corr3p_Proton_EtaDiff_EP_V0A_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0A),WgtEPProton*fWgtCent);
2846  fHist_Corr3p_Proton_EtaDiff_EP_V0C_NN[QAindex][cIndex]->Fill(TMath::Abs(dEta1-dEta2), TMath::Cos(n*dPhi1 + m*dPhi2 - p*PsiNV0C),WgtEPProton*fWgtCent);
2847  }
2848  */
2849  }
2850  //------------------------------------
2851  }
2852 
2853  //cout<<", passes 3 WgtEP = "<<WgtEP<<endl;
2854  }//-------- nested track loop ends ------------------
2855 
2856 
2857 
2858 
2859 
2860 
2861 
2862 
2863 
2864 
2865  //------------- Fill QA histograms ------------
2866  if(TMath::Abs(pVtxZ) < 8.0 && gCharge1 > 0){
2867 
2868  fQAEtaPhiAfterNUA->Fill(dPhi1,dEta1,w1NUA);
2869  /*
2870  if(isPion1){
2871  fQAEtaPhiAfterNUAPion->Fill(dPhi1,dEta1,wNUAPion1);
2872  }
2873  if(isKaon1){
2874  fQAEtaPhiAfterNUAKaon->Fill(dPhi1,dEta1,wNUAKaon1);
2875  }
2876  if(isProton1){
2877  fQAEtaPhiAfterNUAProton->Fill(dPhi1,dEta1,wNUAProton1);
2878  } */
2879  }
2880 
2881 
2882 
2883 
2884 
2885 
2886 
2887 
2888  //-------------- Fill NUA for PID tracks ----------------
2889  if(bFillNUAHistPID) {
2890 
2891  // Pion
2892  if(isPion1){
2893  if(gCharge1>0){
2894  fHist3DEtaPhiVz_Pos_Run[0][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2895  }
2896  else if(gCharge1<0){
2897  fHist3DEtaPhiVz_Neg_Run[0][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2898  }
2899  }
2900 
2901  // Kaon
2902  if(isKaon1){
2903  if(gCharge1>0){
2904  fHist3DEtaPhiVz_Pos_Run[1][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2905  }
2906  else if(gCharge1<0){
2907  fHist3DEtaPhiVz_Neg_Run[1][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2908  }
2909  }
2910 
2911  // proton
2912  if(isProton1){
2913  if(gCharge1>0){
2914  fHist3DEtaPhiVz_Pos_Run[2][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2915  }
2916  else if(gCharge1<0){
2917  fHist3DEtaPhiVz_Neg_Run[2][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2918  }
2919  }
2920 
2921  // Charged
2922  if(gCharge1>0){
2923  fHist3DEtaPhiVz_Pos_Run[3][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2924  }
2925  else if(gCharge1<0){
2926  fHist3DEtaPhiVz_Neg_Run[3][cForNUA]->Fill(pVtxZ,dPhi1,dEta1);
2927  }
2928 
2929  }// Fill NUA for PID or not?
2930 
2931 
2932 
2933 
2934 
2935  //============ PID business starts here =============
2936 
2937  /*
2938  fHistEtaPtAfter->Fill(dEta1,dPt1);
2939 
2940  fHistTPCdEdxvsPAfter->Fill(mom*gCharge1,dEdx1);
2941 
2942  if(TOFmatch>0 && probMis < 0.01){
2943  fHistTOFBetavsPAfter->Fill(mom*gCharge1,beta);
2944  }
2945 
2946  //nSigTOFpion=fPIDResponse->NumberOfSigmasTOF(track, AliPID::kPion);
2947  //nSigTOFkaon=fPIDResponse->NumberOfSigmasTOF(track, AliPID::kKaon);
2948  //nSigTOFproton=fPIDResponse->NumberOfSigmasTOF(track, AliPID::kProton);
2949 
2950  fHistTOFMatchCount->Fill(TOFmatch+2,nSigTOFpion);
2951  fHistTOFMatchCount->Fill(TOFmatch+4,nSigTOFkaon);
2952  fHistTOFMatchCount->Fill(TOFmatch+6,nSigTOFproton);
2953 
2954  if(!TOFmatch || probMis > 0.01){ // I dont want mismatched track in my signal distribution
2955  nSigTOFpion = -99;
2956  nSigTOFkaon = -99;
2957  nSigTOFproton = -99;
2958  }
2959 
2960  //nSigTPCpion = fPIDResponse->NumberOfSigmasTPC(track, AliPID::kPion);
2961  //nSigTPCkaon = fPIDResponse->NumberOfSigmasTPC(track, AliPID::kKaon);
2962  //nSigTPCproton = fPIDResponse->NumberOfSigmasTPC(track, AliPID::kProton);
2963 
2964 
2965  // Fill only 3D nSigma TOF and TPC:
2966  //0=pi, 1=K, 2=Proton
2967  fHistTPCTOFnSigmavsPtAfter[0]->Fill(dPt1,nSigTPCpion,nSigTOFpion);
2968  fHistTPCTOFnSigmavsPtAfter[1]->Fill(dPt1,nSigTPCkaon,nSigTOFkaon);
2969  fHistTPCTOFnSigmavsPtAfter[2]->Fill(dPt1,nSigTPCproton,nSigTOFproton);
2970 
2971 
2972  // Fill nSigTOF for fixed nSigmaTPC Cut
2973  if(TMath::Abs(nSigTPCpion)<=fNSigmaCut){
2974  fHistTOFnSigmavsPtAfter[0]->Fill(dPt1*gCharge1,nSigTOFpion);
2975  fHistPtwithTPCNsigma[0]->Fill(dPt1*gCharge1);
2976  fHistTPCdEdxvsPtPIDAfter[0]->Fill(dPt1,dEdx1);
2977  if(TOFmatch>0 && probMis < 0.01 && beta>0.2){
2978  fHistPtwithTOFSignal[0]->Fill(dPt1*gCharge1);
2979  }
2980  }
2981  if(TMath::Abs(nSigTPCkaon)<=fNSigmaCut){
2982  fHistPtwithTPCNsigma[1]->Fill(dPt1*gCharge1);
2983  fHistTOFnSigmavsPtAfter[1]->Fill(dPt1*gCharge1,nSigTOFkaon);
2984  fHistTPCdEdxvsPtPIDAfter[1]->Fill(dPt1,dEdx1);
2985  if(TOFmatch>0 && probMis < 0.01 && beta>0.2){
2986  fHistPtwithTOFSignal[1]->Fill(dPt1*gCharge1);
2987  }
2988  }
2989  if(TMath::Abs(nSigTPCproton)<=fNSigmaCut){
2990  fHistPtwithTPCNsigma[2]->Fill(dPt1*gCharge1);
2991  fHistTOFnSigmavsPtAfter[2]->Fill(dPt1*gCharge1,nSigTOFproton);
2992  fHistTPCdEdxvsPtPIDAfter[2]->Fill(dPt1,dEdx1);
2993  if(TOFmatch>0 && probMis < 0.01 && beta>0.2){
2994  fHistPtwithTOFSignal[2]->Fill(dPt1*gCharge1);
2995  }
2996  }
2997 
2998 
2999  // Fill dEdx Histograms
3000  if(TMath::Abs(nSigTPCpion)<=fNSigmaCut){
3001  fHistPtwithTPCNsigma[0]->Fill(dPt1*gCharge1);
3002  fHistTPCdEdxvsPtPIDAfter[0]->Fill(dPt1,dEdx1);
3003  }
3004  if(TMath::Abs(nSigTPCkaon)<=fNSigmaCut){
3005  fHistPtwithTPCNsigma[1]->Fill(dPt1*gCharge1);
3006  fHistTPCdEdxvsPtPIDAfter[1]->Fill(dPt1,dEdx1);
3007  }
3008  if(TMath::Abs(nSigTPCproton)<=fNSigmaCut){
3009  fHistPtwithTPCNsigma[2]->Fill(dPt1*gCharge1);
3010  fHistTPCdEdxvsPtPIDAfter[2]->Fill(dPt1,dEdx1);
3011  }
3012 
3013 
3014 
3015  //========> nSigmaTOF distribution for circular cut <===============
3016  //if(TMath::Sqrt(nSigTPCpion*nSigTPCpion+nSigTOFpion*nSigTOFpion)<=fNSigmaCut){
3017  //fHistTOFnSigmavsPtAfter[0]->Fill(dPt1*gCharge1,nSigTOFpion);
3018  //if(TOFmatch>0 && probMis < 0.01 && beta>0.2){
3019  //fHistPtwithTOFSignal[0]->Fill(dPt1*gCharge1);
3020  //}
3021  //}
3022  //if(TMath::Sqrt(nSigTPCkaon*nSigTPCkaon+nSigTOFkaon*nSigTOFkaon)<=fNSigmaCut){
3023  //fHistTOFnSigmavsPtAfter[1]->Fill(dPt1*gCharge1,nSigTOFkaon);
3024  //if(TOFmatch>0 && probMis < 0.01 && beta>0.2){
3025  //fHistPtwithTOFSignal[1]->Fill(dPt1*gCharge1);
3026  //}
3027  //}
3028  //if(TMath::Sqrt(nSigTPCproton*nSigTPCproton+nSigTOFproton*nSigTOFproton)<=fNSigmaCut){
3029  //fHistTOFnSigmavsPtAfter[2]->Fill(dPt1*gCharge1,nSigTOFproton);
3030  //if(TOFmatch>0 && probMis < 0.01 && beta>0.2){
3031  //fHistPtwithTOFSignal[2]->Fill(dPt1*gCharge1);
3032  //}
3033  //}
3034 
3035 
3036  //------------------------------------------------------------
3037 
3038  // TOF matching is not needed from data as done by Davide using MC.
3039 
3040  if(!TOFmatch || probMis > 0.01 || beta>0.2) continue;
3041 
3042  // nSigmaTPC distribution for Fixed nSigmaTOF cut
3043  if(TMath::Abs(nSigTOFpion)<=fNSigmaCut){
3044  fHistTPCnSigmavsPtAfter[0]->Fill(dPt1*gCharge1,nSigTPCpion);
3045  fHistPtwithTOFmasscut[0]->Fill(dPt1*gCharge1);
3046  }
3047  if(TMath::Abs(nSigTOFkaon)<=fNSigmaCut){
3048  fHistTPCnSigmavsPtAfter[1]->Fill(dPt1*gCharge1,nSigTPCkaon);
3049  fHistPtwithTOFmasscut[1]->Fill(dPt1*gCharge1);
3050  }
3051  if(TMath::Abs(nSigTOFproton)<=fNSigmaCut){
3052  fHistTPCnSigmavsPtAfter[2]->Fill(dPt1*gCharge1,nSigTPCproton);
3053  fHistPtwithTOFmasscut[2]->Fill(dPt1*gCharge1);
3054  } */
3055 
3056 
3057 
3058  }
3059 //===================== track loop ends ============================
3060 
3061 
3062 
3063 
3064 
3065 
3066 
3067  //---------- TPC Event Planes -------
3068  //Sub events //
3069  sumTPCQn2x[0] = sumTPCQn2x[0]/SumWEtaNeg;
3070  sumTPCQn2y[0] = sumTPCQn2y[0]/SumWEtaNeg;
3071  sumTPCQn2x[1] = sumTPCQn2x[1]/SumWEtaPos;
3072  sumTPCQn2y[1] = sumTPCQn2y[1]/SumWEtaPos;
3073 
3074  sumTPCQn2x[2] = sumTPCQn2x[2]/SumWEtaFull;
3075  sumTPCQn2y[2] = sumTPCQn2y[2]/SumWEtaFull;
3076 
3077  /*
3078  sumTPCQn3x[0] = sumTPCQn3x[0]/SumWEtaNeg;
3079  sumTPCQn3y[0] = sumTPCQn3y[0]/SumWEtaNeg;
3080  sumTPCQn4x[0] = sumTPCQn4x[0]/SumWEtaNeg;
3081  sumTPCQn4y[0] = sumTPCQn4y[0]/SumWEtaNeg;
3082 
3083  sumTPCQn3x[1] = sumTPCQn3x[1]/SumWEtaPos;
3084  sumTPCQn3y[1] = sumTPCQn3y[1]/SumWEtaPos;
3085  sumTPCQn4x[1] = sumTPCQn4x[1]/SumWEtaPos;
3086  sumTPCQn4y[1] = sumTPCQn4y[1]/SumWEtaPos;
3087 
3088  sumTPCQn3x[2] = sumTPCQn3x[2]/SumWEtaFull;
3089  sumTPCQn3y[2] = sumTPCQn3y[2]/SumWEtaFull;
3090  sumTPCQn4x[2] = sumTPCQn4x[2]/SumWEtaFull;
3091  sumTPCQn4y[2] = sumTPCQn4y[2]/SumWEtaFull;
3092  */
3093 
3094  if(gPsiN==3){
3095  //Enable periodicity PsiN directly:
3096  PsiNTPCC = (1.0/gPsiN)*( TMath::ATan2(sumTPCQn2y[0],sumTPCQn2x[0]) + TMath::Pi() ) ; // negetive eta
3097  PsiNTPCA = (1.0/gPsiN)*( TMath::ATan2(sumTPCQn2y[1],sumTPCQn2x[1]) + TMath::Pi() ) ; // positive eta
3098  PsiNTPCF = (1.0/gPsiN)*( TMath::ATan2(sumTPCQn2y[2],sumTPCQn2x[2]) + TMath::Pi() ) ; // FUll TPC eta
3099  }
3100  else{
3101  PsiNTPCC = (1.0/gPsiN)*( TMath::ATan2(sumTPCQn2y[0],sumTPCQn2x[0]) ) ; // negetive eta
3102  if(PsiNTPCC<0.) PsiNTPCC += 2*TMath::Pi()/gPsiN;
3103  PsiNTPCA = (1.0/gPsiN)*( TMath::ATan2(sumTPCQn2y[1],sumTPCQn2x[1]) ) ; // positive eta
3104  if(PsiNTPCA<0.) PsiNTPCA += 2*TMath::Pi()/gPsiN;
3105  PsiNTPCF = (1.0/gPsiN)*( TMath::ATan2(sumTPCQn2y[2],sumTPCQn2x[2]) ) ; // FUll TPC eta
3106  if(PsiNTPCF<0.) PsiNTPCF += 2*TMath::Pi()/gPsiN;
3107  }
3108 
3109  fHTPCAEventPlaneVsCent->Fill(EvtCent,PsiNTPCA);
3110  fHTPCCEventPlaneVsCent->Fill(EvtCent,PsiNTPCC);
3111  fHTPCEventPlaneVsCent->Fill(EvtCent, PsiNTPCF);
3112 
3113 
3114 
3115 
3116 
3117  //-------------- vs Centrality -----------------
3118  //V0A-V0C
3119  fHist_Reso2n_EP_Norm_Det[QAindex][0]->Fill(EvtCent, TMath::Cos(gPsiN*(PsiNV0A-PsiNV0C)), fWgtCent);
3120  //V0A-TPC
3121  fHist_Reso2n_EP_Norm_Det[QAindex][1]->Fill(EvtCent, TMath::Cos(gPsiN*(PsiNV0A-PsiNTPCF)), fWgtCent);
3122  //V0C-TPC
3123  fHist_Reso2n_EP_Norm_Det[QAindex][2]->Fill(EvtCent, TMath::Cos(gPsiN*(PsiNV0C-PsiNTPCF)), fWgtCent);
3124  //TPCa -TPCc
3125  fHist_Reso2n_EP_Norm_Det[QAindex][3]->Fill(EvtCent, TMath::Cos(gPsiN*(PsiNTPCA-PsiNTPCC)), fWgtCent);
3126 
3127 
3128  //--------- Store TPC-Qn for Recenter ---------
3129  fTPCAQ2xVsCentRun->Fill(EvtCent,sumTPCQn2x[0]);
3130  fTPCAQ2yVsCentRun->Fill(EvtCent,sumTPCQn2y[0]);
3131  fTPCCQ2xVsCentRun->Fill(EvtCent,sumTPCQn2x[1]);
3132  fTPCCQ2yVsCentRun->Fill(EvtCent,sumTPCQn2y[1]);
3133 
3134  fTPCFQ2xVsCentRun->Fill(EvtCent,sumTPCQn2x[2]);
3135  fTPCFQ2yVsCentRun->Fill(EvtCent,sumTPCQn2y[2]);
3136  /*
3137  fTPCAQ3xVsCentRun->Fill(EvtCent,sumTPCQn3x[0]);
3138  fTPCAQ3yVsCentRun->Fill(EvtCent,sumTPCQn3y[0]);
3139  fTPCCQ3xVsCentRun->Fill(EvtCent,sumTPCQn3x[1]);
3140  fTPCCQ3yVsCentRun->Fill(EvtCent,sumTPCQn3y[1]);
3141 
3142  fTPCAQ4xVsCentRun->Fill(EvtCent,sumTPCQn4x[0]);
3143  fTPCAQ4yVsCentRun->Fill(EvtCent,sumTPCQn4y[0]);
3144  fTPCCQ4xVsCentRun->Fill(EvtCent,sumTPCQn4x[1]);
3145  fTPCCQ4yVsCentRun->Fill(EvtCent,sumTPCQn4y[1]);
3146  */
3147 
3148 
3149 
3150 
3151 
3152 
3153 
3154 
3155 
3156 
3157 
3158  PostData(1,fListHist);
3159 
3160  fHistEventCount->Fill(14.5); //15th bin is last one
3161  stepCount++;
3162 
3163 
3164 
3165 
3166 
3167 
3168  //if(multPOI2nd!=(multPOI1st-1))
3169  //cout<<"mismatched "<< "\tPOIs1st = "<<multPOI1st<<"\tPOIs2nd = "<<multPOI2nd <<" for Event = "<<fEventCount<<endl;
3170 
3171  //if(multPOI1st!=multEtaFull)
3172  //cout<<"\n Mismatched track "<< "\t multEtaFull = "<<multEtaFull<<"\tPOIs1st = "<<multPOI1st<<" for Event = "<<fEventCount<<endl;
3173 
3174  //if(fEventCount%10==0)
3175  //cout<<"Ev = "<<fEventCount<<"\tMult = "<<multEtaFull<<"\tPOIs1st = "<<multPOI1st<<"\t POI2nd = "<< multPOI2nd <<endl;
3176 
3177  //if(multEtaFull>500)
3178  //cout<<"Ev = "<<fEventCount<<" ntracks = "<<ntracks<<"\tMult = "<<multEtaFull<<"\t RealTime = "<<watch.RealTime()<<"\t CPUTime = "<< watch.CpuTime() <<endl;
3179  //cout<<"Ev = "<<fEventCount<<" cent = "<<EvtCent<<"\tWEtaNeg = "<<SumWEtaNeg<<"\tnegMult = "<<multEtaNeg<<"\tWEtaPos = "<<SumWEtaPos<<"\tposMult = "<<multEtaPos<<endl;
3180  //watch.Stop();
3181 
3182 }//================ UserExec ==============
3183 
3184 
3185 
3186 
3187 
3188 
3189 
3190 
3191 
3192 
3193 
3194 
3195 
3196 
3197 
3198 
3199 
3200 
3201 
3202 
3203 
3204 
3206 
3207 
3208 
3209 
3210 double AliAnalysisTaskCMEV0PID::GetWDist(const AliVVertex* v0, const AliVVertex* v1)
3211 {
3212  // calculate sqrt of weighted distance to other vertex
3213  if (!v0 || !v1) {
3214  AliDebug(2,"\n\n ::GetWDist => One of vertices is not valid\n\n");
3215  return 0;
3216  }
3217  static TMatrixDSym vVb(3);
3218  double dist = -1;
3219  double dx = v0->GetX()-v1->GetX();
3220  double dy = v0->GetY()-v1->GetY();
3221  double dz = v0->GetZ()-v1->GetZ();
3222  double cov0[6],cov1[6];
3223  v0->GetCovarianceMatrix(cov0);
3224  v1->GetCovarianceMatrix(cov1);
3225  vVb(0,0) = cov0[0]+cov1[0];
3226  vVb(1,1) = cov0[2]+cov1[2];
3227  vVb(2,2) = cov0[5]+cov1[5];
3228  vVb(1,0) = vVb(0,1) = cov0[1]+cov1[1];
3229  vVb(0,2) = vVb(1,2) = vVb(2,0) = vVb(2,1) = 0.;
3230  vVb.InvertFast();
3231  if (!vVb.IsValid()) {
3232  AliDebug(2,"Singular Matrix\n");
3233  return dist;
3234  }
3235  dist = vVb(0,0)*dx*dx + vVb(1,1)*dy*dy + vVb(2,2)*dz*dz
3236  + 2*vVb(0,1)*dx*dy + 2*vVb(0,2)*dx*dz + 2*vVb(1,2)*dy*dz;
3237  return dist>0 ? TMath::Sqrt(dist) : -1;
3238 }
3239 
3240 
3242  { // check for multi-vertexer pile-up
3243  const int kMinPlpContrib = 5;
3244  const double kMaxPlpChi2 = 5.0;
3245  const double kMinWDist = 15;
3246 
3247  const AliVVertex* vtPrm = 0;
3248  const AliVVertex* vtPlp = 0;
3249 
3250  int nPlp = 0;
3251 
3252  if(!(nPlp=faod->GetNumberOfPileupVerticesTracks()))
3253  return kFALSE;
3254 
3255  vtPrm = faod->GetPrimaryVertex();
3256  if(vtPrm == faod->GetPrimaryVertexSPD())
3257  return kTRUE; // there are pile-up vertices but no primary
3258 
3259  //int bcPrim = vtPrm->GetBC();
3260 
3261  for(int ipl=0;ipl<nPlp;ipl++) {
3262  vtPlp = (const AliVVertex*)faod->GetPileupVertexTracks(ipl);
3263  if (vtPlp->GetNContributors() < kMinPlpContrib) continue;
3264  if (vtPlp->GetChi2perNDF() > kMaxPlpChi2) continue;
3265  //int bcPlp = vtPlp->GetBC();
3266  //if (bcPlp!=AliVTrack::kTOFBCNA && TMath::Abs(bcPlp-bcPrim)>2)
3267  // return kTRUE; // pile-up from other BC
3268 
3269  double wDst = GetWDist(vtPrm,vtPlp);
3270  if (wDst<kMinWDist) continue;
3271 
3272  return kTRUE; // pile-up: well separated vertices
3273  }
3274  return kFALSE;
3275 }
3276 
3277 
3278 
3280  /*if(bApplyMCcorr){
3281  if(!gGrid){
3282  TGrid::Connect("alien://");
3283  }
3284  if(!mfileFBHijing){
3285  mfileFBHijing = TFile::Open(sMCfilePath,"READ");
3286  fListFBHijing = dynamic_cast<TList*>(mfileFBHijing->FindObjectAny("fMcEffiHij")); */
3287 
3288  if(fListFBHijing) {
3289  cout<<"\n =========> Info: Using MC efficiency correction Map <=========== "<<endl;
3290  for(int i=0;i<10;i++) {
3291  fFB_Efficiency_Cent[i] = (TH1D *) fListFBHijing->FindObject(Form("eff_unbiased_%d",i));
3292  }
3293  //PID :
3294  /*
3295  for(int i=0;i<10;i++) {
3296  fFB_Efficiency_Pion_Cent[i] = (TH1D *) fListFBHijing->FindObject(Form("eff_unbiased_Pion_%d",i));
3297  fFB_Efficiency_Kaon_Cent[i] = (TH1D *) fListFBHijing->FindObject(Form("eff_unbiased_Kaon_%d",i));
3298  fFB_Efficiency_Proton_Pos_Cent[i] = (TH1D *) fListFBHijing->FindObject(Form("eff_unbiased_ProtonPos_%d",i));
3299  fFB_Efficiency_Proton_Neg_Cent[i] = (TH1D *) fListFBHijing->FindObject(Form("eff_unbiased_ProtonNeg_%d",i));
3300  }
3301  */
3302  //------- Fill the flags: --------------
3304  fHistTaskConfigParameters->SetBinContent(19,1);
3305  }
3306  /*
3307  if(fFB_Efficiency_Pion_Cent[0] && fFB_Efficiency_Pion_Cent[4]){
3308  fHistTaskConfigParameters->SetBinContent(20,1);
3309  }
3310  if(fFB_Efficiency_Kaon_Cent[0] && fFB_Efficiency_Kaon_Cent[4]){
3311  fHistTaskConfigParameters->SetBinContent(21,1);
3312  }
3313  if(fFB_Efficiency_Proton_Pos_Cent[0] && fFB_Efficiency_Proton_Neg_Cent[0]){
3314  fHistTaskConfigParameters->SetBinContent(22,1);
3315  } */
3316  }
3317  else if(!fListFBHijing){
3318  std::cout<<"\n\n !!!!**** Warning: FB Efficiency File/List not found Use Weight = 1.0 *****\n\n"<<std::endl;
3319  //exit(1);
3320  }
3321 }
3322 
3323 
3324 
3325 
3327 
3328  Int_t cIndex = 0;
3329 
3330  if(fCent<5.0) {
3331  cIndex = 0;
3332  }
3333  else if(fCent>=5.0 && fCent<10){
3334  cIndex = 1;
3335  }
3336  else if(fCent>=10.0) {
3337  cIndex = abs(fCent/10.0)+1;
3338  }
3339  return cIndex;
3340 }
3341 
3343  //std::cout<<" centrality outlier function called "<<std::endl;
3344  Float_t fMeanTPC[100] = {2902.95,2758.33,2642.78,2536.67,2435.37,2340.06,2248.44,2163.71,2080.49,2001.54,1925.86,1852.64,1781.97,1715.56,1650.53,1587.23,1527.51,1468.19,1412.73,1357.86,1305.35,1254.33,1205.57,1157.28,1111.53,1066.42,1023.15,981.594,940.795,901.766,863.651,826.183,790.53,756.358,722.654,690.513,659.443,628.807,599.748,571.664,544.446,518.042,492.369,468.072,444.694,422.487,400.104,379.129,359.147,339.62,320.817,302.788,285.791,269.015,253.688,238.671,224.039,209.932,196.915,184.647,172.76,161.381,150.395,140.288,131.033,121.58,113.112,104.938,97.3078,90.2178,83.5974,77.2645,70.7126,65.4424,60.1404,55.5644,50.8314,46.3761,43.024,38.625,35.3435,32.2304,29.4192,26.821,24.3303,21.9332,19.4215,16.7163,14.9414,13.1092,0.};
3345 
3346  Float_t fSigmaTPC[100] = {122.209,107.901,103.452,100.498,97.7403,94.7845,93.2543,90.0548,88.1106,85.7382,84.0812,82.2978,80.3817,78.6002,77.3448,75.5086,73.6842,71.9733,70.3447,69.1999,67.878,66.3511,65.0406,63.4866,62.4409,60.7899,59.1328,58.426,56.8618,55.8871,54.1031,53.4959,52.0482,51.0441,49.6218,48.7646,47.5166,46.5247,45.0727,44.4311,43.4531,42.0404,41.0238,40.1384,39.2588,38.2461,36.5951,36.0552,35.3727,33.7883,32.7167,32.4486,31.3709,30.3444,29.505,28.5139,27.4471,26.5359,25.9506,25.127,24.3797,23.2985,22.279,21.4698,20.781,20.8193,19.9509,18.8036,17.9145,16.961,16.7375,15.852,14.9324,14.7663,13.5969,13.4533,12.3067,12.7835,11.7283,10.6758,10.6676,10.6492,9.04614,8.89065,8.66093,8.50997,7.98812,6.91087,7.12045,7.29593,0.};
3347 
3348  for(int i=0;i<90;i++) {
3349  hCentvsTPCmultCuts->SetBinContent(i+1,1,fMeanTPC[i]);
3350  hCentvsTPCmultCuts->SetBinContent(i+1,2,fSigmaTPC[i]);
3351  }
3352 }
3353 
3354 
3355 
3356 
3357 
3358 
3360 
3361  fHistTaskConfigParameters = new TH1F("fHistTaskConfigParameters","Task parameters",25,0,25);
3362  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(1,"FilterBit");
3363  fHistTaskConfigParameters->SetBinContent(1,fFilterBit);
3364  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(2,"n#sigmaTPC");
3365  fHistTaskConfigParameters->SetBinContent(2,fNSigmaCut);
3366  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(3,"MinPt");
3367  fHistTaskConfigParameters->SetBinContent(3,fMinPtCut);
3368  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(4,"MaxPt");
3369  fHistTaskConfigParameters->SetBinContent(4,fMaxPtCut);
3370  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(5,"MinEta");
3371  fHistTaskConfigParameters->SetBinContent(5,fMinEtaCut);
3372  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(6,"MaxEta");
3373  fHistTaskConfigParameters->SetBinContent(6,fMaxEtaCut);
3374 
3375  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(11,"CentralityMin");
3377  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(12,"CentralityMax");
3379 
3380  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(13,"VertexMin(cm)");
3381  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(14,"VertexMax(cm)");
3382 
3383 
3384  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(15,"NUA Charge");
3385  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(16,"NUA Pion");
3386  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(17,"NUA Kaon");
3387  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(18,"NUA Proton");
3388 
3389  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(19,"MC Charge");
3390  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(20,"MC Pion");
3391  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(21,"MC Kaon");
3392  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(22,"MC Proton");
3393 
3394  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(23,"V0 Gain Eq.");
3395  fHistTaskConfigParameters->GetXaxis()->SetBinLabel(24,"V0 Qn Recenter");
3396 
3397 
3399 
3400 
3401 
3402  fHistPileUpCount = new TH1F("fHistPileUpCount", "fHistPileUpCount", 15, 0., 15.);
3403  fHistPileUpCount->GetXaxis()->SetBinLabel(1,"plpMV");
3404  fHistPileUpCount->GetXaxis()->SetBinLabel(2,"fromSPD");
3405  fHistPileUpCount->GetXaxis()->SetBinLabel(3,"RefMultComb08");
3406  fHistPileUpCount->GetXaxis()->SetBinLabel(4,"IncompleteDAQ");
3407  fHistPileUpCount->GetXaxis()->SetBinLabel(5,"abs(V0M-CL1)>5.0");
3408  fHistPileUpCount->GetXaxis()->SetBinLabel(6,"missingVtx");
3409  fHistPileUpCount->GetXaxis()->SetBinLabel(7,"inconsistentVtx");
3410  Int_t puConst = fPileUpConstParm;
3411  fHistPileUpCount->GetXaxis()->SetBinLabel(8,Form("multESDTPCDif>%d",puConst));
3412  fHistPileUpCount->GetXaxis()->SetBinLabel(9,Form("multGlobTPCDif>%d",puConst));
3413  fHistPileUpCount->GetXaxis()->SetBinLabel(10,"PileUpMultSelTask");
3415 
3416 
3417  fHistMultSelPUCount = new TH1F("fHistMultSelPileUpCount", "no of PU Event from MultSelTask", 10, 0., 10);
3418  fHistMultSelPUCount->GetXaxis()->SetBinLabel(1,"PileUp");
3419  fHistMultSelPUCount->GetXaxis()->SetBinLabel(2,"PileUpMV");
3420  fHistMultSelPUCount->GetXaxis()->SetBinLabel(3,"PileUpMultBins");
3421  fHistMultSelPUCount->GetXaxis()->SetBinLabel(4,"InconsistentVtx");
3422  fHistMultSelPUCount->GetXaxis()->SetBinLabel(5,"TrackletVsCluster");
3423  fHistMultSelPUCount->GetXaxis()->SetBinLabel(6,"IncompleteDAQ");
3424  fHistMultSelPUCount->GetXaxis()->SetBinLabel(7,"NotGoodVertex2016");
3426 
3427 
3428  fHistEventCount = new TH1F("fHistEventCount","Event counts",15,0,15);
3429  fHistEventCount->GetXaxis()->SetBinLabel(1,"Called UserExec()");
3430  fHistEventCount->GetXaxis()->SetBinLabel(2,"Called Exec()");
3431  fHistEventCount->GetXaxis()->SetBinLabel(3,"AOD Exist");
3432  fHistEventCount->GetXaxis()->SetBinLabel(4,"Vz < 10");
3433  fHistEventCount->GetXaxis()->SetBinLabel(5,Form("%2.0f<Cent<%2.0f",fCentralityPercentMin,fCentralityPercentMax));
3434  fHistEventCount->GetXaxis()->SetBinLabel(6,"noAODtrack > 2 ");
3435  fHistEventCount->GetXaxis()->SetBinLabel(7,"TPC vs Global");
3436  fHistEventCount->GetXaxis()->SetBinLabel(8,"TPC128 vs ESD");
3437  fHistEventCount->GetXaxis()->SetBinLabel(9,"Cent vs TPC");
3438  fHistEventCount->GetXaxis()->SetBinLabel(10,"mult eta+/- > 2");
3439  fHistEventCount->GetXaxis()->SetBinLabel(11,"centCL1 < 90");
3440  fHistEventCount->GetXaxis()->SetBinLabel(15,"Survived Events");
3441  fListHist->Add(fHistEventCount);
3442 
3443  //fHistEventCount->Fill(1);
3444 
3445 }
3446 
3447 
3448 
3449 
3450 
3452 {
3453  /*
3454  TList *fListAppliedV0Corr = new TList();
3455  fListAppliedV0Corr->SetName("fListAppliedV0Corr");
3456  fListAppliedV0Corr->SetOwner(kTRUE);
3457  fListHist->Add(fListAppliedV0Corr);*/ // to be Tested later.
3458 
3459  if(fListV0MCorr){
3460  fHCorrectV0M = (TH1D *) fListV0MCorr->FindObject(Form("fHistV0Gain_Run%d",run));
3461  fHAvgerageQnV0A = (TH2D *) fListV0MCorr->FindObject(Form("fHistAvgQnV0A_Run%d",run));
3462  fHAvgerageQnV0C = (TH2D *) fListV0MCorr->FindObject(Form("fHistAvgQnV0C_Run%d",run));
3463 
3464  //Now the wgts for centrality is added inside VOM gain file.
3465  fHCentWeightForRun = (TH1D *) fListV0MCorr->FindObject(Form("fHistCentWeight_Run%d",run));
3466 
3467  if(fHCentWeightForRun){
3468  cout<<"\n =========== Info:: Using Centrality wgts for run = "<<run<<"============"<<endl;
3469  }
3470  else{
3471  cout<<"\n =========== Info:: No Centrality wgt. correction.!! for run = "<<run<<"============"<<endl;
3472  }
3473 
3474  if(fHCorrectV0M){
3475  cout<<"\n =========== Info:: Setting up V0 gain correction for run = "<<run<<"============"<<endl;
3476  fHistTaskConfigParameters->SetBinContent(23,1);
3477  }
3478  else{
3479  cout<<"\n =========== Info:: No V0 gain correction..!!! for run = "<<run<<"============"<<endl;
3480  }
3482  cout<<"\n =========== Info:: Setting up V0 <Qn> recentering for run = "<<run<<"============"<<endl;
3483  fHistTaskConfigParameters->SetBinContent(24,1);
3484  }
3485  else{
3486  cout<<"\n =========== Info:: No V0 <Qn> recentering..!!! for run = "<<run<<"============"<<endl;
3487  }
3488  }
3489  else{
3490  cout<<"\n ======== !!!! Error:: List For V0-Gain and Qn not found for run "<<run<<"============"<<endl;
3491  }
3492 
3493 }
3494 
3495 
3496 
3497 
3498 
3499 
3500 
3502 {
3503 
3504  if(fListNUACorr){
3505  for(int i=0;i<5;i++){
3506  fHCorrectNUApos[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Pos_Cent%d_Run%d",i,run));
3507  fHCorrectNUAneg[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Neg_Cent%d_Run%d",i,run));
3508  }
3509  if(fHCorrectNUApos[3] && fHCorrectNUAneg[3]){
3510  cout<<"\n=========== Info:: Setting up NUA corrections for run "<<run<<"============"<<endl;
3511  fHistTaskConfigParameters->SetBinContent(15,1);
3512  }
3513  }
3514  else {
3515  printf("\n ******** Warning: No NUA Correction for Charge Particle in run %d, Use Wgt = 1.0 ********* \n",run);
3516  }
3517 
3518  //=================== PID: ==========================
3519  if(fListNUACorr){
3520  /*
3521  for(int i=0;i<5;i++){
3522  fHCorrectNUAposPion[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Pion_Pos_Cent%d_Run%d",i,run)); //
3523  fHCorrectNUAnegPion[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Pion_Neg_Cent%d_Run%d",i,run));
3524  if(fHCorrectNUAposPion[3] && fHCorrectNUAnegPion[3]) fHistTaskConfigParameters->SetBinContent(16,1);
3525 
3526  fHCorrectNUAposKaon[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Kaon_Pos_Cent%d_Run%d",i,run)); //
3527  fHCorrectNUAnegKaon[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Kaon_Neg_Cent%d_Run%d",i,run));
3528  if(fHCorrectNUAposKaon[3] && fHCorrectNUAnegKaon[3]) fHistTaskConfigParameters->SetBinContent(17,1);
3529 
3530  fHCorrectNUAposProton[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Proton_Pos_Cent%d_Run%d",i,run));
3531  fHCorrectNUAnegProton[i] = (TH3F *) fListNUACorr->FindObject(Form("fHist_NUA_VzPhiEta_Proton_Neg_Cent%d_Run%d",i,run));
3532  if(fHCorrectNUAposProton[3] && fHCorrectNUAnegProton[3]) fHistTaskConfigParameters->SetBinContent(18,1);
3533  }
3534  if(fHCorrectNUAposPion[3] && fHCorrectNUAposKaon[3] && fHCorrectNUAposProton[3]) {
3535  cout<<"\n=========== Info:: Setting up --> PID NUA corrections for run = "<<run<<"============"<<endl;
3536  }
3537  else{
3538  cout<<"\n=========== WARNING :: PID NUA corrections NOT found for run = "<<run<<"============"<<endl;
3539  }*/
3540  }
3541  else {
3542  printf("\n ******** Error/Warning: No NUA Correction found for PID for run %d, Use PID Wgt = 1.0 ********* \n",run);
3543  }
3544 
3545 }
3546 
TProfile * fHist_Corr3p_pTSum_EP_V0C_PP[2][6]
Int_t GetCentralityScaled0to10(Double_t fCent)
double Double_t
Definition: External.C:58
Definition: External.C:260
TH2D * fHAvgerageQnV0A
To read Gain Correction file.
TProfile * fHist_Corr3p_pTDiff_EP_V0C_PP[2][6]
Definition: External.C:236
const char * title
Definition: MakeQAPdf.C:27
TProfile * fHist_Corr3p_pTSum_EP_V0A_PN[2][6]
TProfile * fHist_Corr3p_pTSum_EP_V0C_NN[2][6]
AliMultSelection * fMultSelection
PID response Handler.
centrality
char Char_t
Definition: External.C:18
TProfile * fHist_Corr3p_EtaDiff_EP_V0A_PN[2][6]
TH1F * fHistPtwithTPCNsigma[3]
last in the list
TProfile * fHist_Reso2n_EP_Norm_Det[2][4]
TH3F * fHCorrectNUAnegKaon[5]
5 centrality bin, read NUA from file
TH2F * fHistEtaPtAfter
Eta-Pt acceptance.
TProfile * fHist_Corr3p_EP_Norm_NN[2][4]
TProfile * fHist_Corr3p_pTDiff_EP_V0C_PN[2][6]
TProfile * fHist_Corr3p_EtaDiff_EP_V0C_PP[2][6]
TH3F * fHCorrectNUAposPion[5]
5 centrality bin, read NUA from file
TH1F * fCentDistAfter
without PileUp cut
TProfile * fHist_Corr3p_EtaDiff_EP_V0A_NN[2][6]
TProfile * fHist_Corr3p_pTSum_EP_V0C_PN[2][6]
TH1D * fHCentWeightForRun
V0C Average <Qn>, n=2,3.
TProfile * fHist_Corr3p_EtaDiff_EP_V0A_PP[2][6]
int Int_t
Definition: External.C:63
float Float_t
Definition: External.C:68
TH3F * fHist3DEtaPhiVz_Neg_Run[4][5]
4 particle 5 centrality bin
Definition: External.C:228
Definition: External.C:212
TH3F * fHCorrectNUAnegProton[5]
5 centrality bin, read NUA from file
TH3F * fHCorrectNUAposKaon[5]
5 centrality bin, read NUA from file
TH1D * fHCorrectV0M
with PileUp cut
TH3F * fHCorrectNUAposProton[5]
5 centrality bin, read NUA from file
TH2F * fHistTPCvsGlobalMultBefore
Eta-Pt acceptance.
TProfile * fHist_Corr3p_pTDiff_EP_V0A_PP[2][6]
TH1F * fCentDistBefore
To fill VOM multiplicity.
TList * fListHist
Event selection.
virtual void UserExec(Option_t *)
TH1F * fHistPileUpCount
Task input parameters FB / cut values etc.
Bool_t PileUpMultiVertex(const AliAODEvent *faod)
double GetWDist(const AliVVertex *v0, const AliVVertex *v1)
TH1D * fFB_Efficiency_Cent[10]
4 particle 5 centrality bin
TProfile * fHist_Corr3p_EP_Norm_PN[2][4]
5 centrality bin, read NUA from file
TProfile * fHist_Corr3p_EtaDiff_EP_V0C_PN[2][6]
TProfile * fHist_Corr3p_pTSum_EP_V0A_NN[2][6]
TProfile * fHist_Corr3p_pTDiff_EP_V0C_NN[2][6]
TH2F * fQAEtaPhiAfterNUA
Event weights for non-flat centrality.
const char Option_t
Definition: External.C:48
TH3F * fHCorrectNUAneg[5]
5 centrality bin, read NUA from file
TH3F * fHCorrectNUAnegPion[5]
5 centrality bin, read NUA from file
bool Bool_t
Definition: External.C:53
TProfile * fHist_Corr3p_pTDiff_EP_V0A_PN[2][6]
TProfile * fHist_Corr3p_pTDiff_EP_V0A_NN[2][6]
TH2F * fV0MultChVsRun
Full Event plane.
TProfile * fHist_Corr3p_pTSum_EP_V0A_PP[2][6]
TProfile * fHist_Corr3p_EtaDiff_EP_V0C_NN[2][6]
TH2D * fHAvgerageQnV0C
V0A Average <Qn>, n=2,3.
TProfile * fHist_Corr3p_EP_Norm_PP[2][4]