]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGCF/FLOW/CME/AliAnalysisTaskCMEv2A.cxx
changes to be able to run on train:
[u/mrichter/AliRoot.git] / PWGCF / FLOW / CME / AliAnalysisTaskCMEv2A.cxx
CommitLineData
545d70ee 1#include <iostream>
2#include <cstdlib>
3#include <sys/time.h>
4
5// ROOT classes
6#include "TChain.h"
7#include "TTree.h"
8#include "TList.h"
9#include "TH1F.h"
10#include "TH2F.h"
11#include "TProfile.h"
12#include "TExMap.h"
13#include "TRandom3.h"
14#include "TMath.h"
15
16// Alice analysis base class
17#include "AliAnalysisTaskSE.h"
18
19// Alice analysis additional classes
20#include "AliAnalysisManager.h"
21//#include "AliInputEventHandler.h"
22
23// Alice AOD classes
24#include "AliAODInputHandler.h"
25//#include "AliAODHandler.h"
26#include "AliAODEvent.h"
27#include "AliAODVertex.h"
28#include "AliAODVZERO.h"
29
30// Alice classes
31#include "AliCentrality.h"
32#include "AliEventplane.h"
33#include "AliAnalysisUtils.h"
34
35// Alice MC classes
36#include "AliMCEvent.h"
37//#include "AliMCEventHandler.h"
38#include "AliAODMCParticle.h"
39
40// Alice "V" classes
41#include "AliVParticle.h"
42//#include "AliVEvent.h"
43//#include "AliVVertex.h"
44//#include "AliVVZERO.h"
45
46// Alice PID classes
47//#include "AliAODPid.h"
48//#include "AliAODpidUtil.h"
49#include "AliPID.h"
50//#include "AliPIDCombined.h"
51#include "AliPIDResponse.h"
52
53
54// This class
55#include "AliAnalysisTaskCMEv2A.h"
56
57ClassImp(AliAnalysisTaskCMEv2A); // import class inheriting from TObject
58
59using std::cout;
60using std::endl;
61
62const float pi = 3.141592653589793;
63
64// function prototype...
65float GetDPhiStar(float phi1, float pt1, float charge1, float phi2, float pt2, float charge2, float radius, float bSign);
66
67
68//------------------------------------------------------------------
69AliAnalysisTaskCMEv2A::AliAnalysisTaskCMEv2A() : AliAnalysisTaskSE()
70{
71 // Default class constructor
72
73 cout<<"Default class constructor called. Prepare for fun analysis time!"<<endl;
74 this->SetParameters();
75
76}
77
78
79
80//--------------------------------------------------------------------------------------
81AliAnalysisTaskCMEv2A::AliAnalysisTaskCMEv2A(const char *name) : AliAnalysisTaskSE(name)
82{
83 // Class Constructor with name
84
85 cout<<"User defined class constructor called. Prepare for fun analysis time!"<<endl;
86 this->SetParameters();
87
88 // Define input and output slots here
89 // Input slot #0 works with a TChain
90 DefineInput(0,TChain::Class());
91
92 // Output slot #0 is reserved by the base class for AOD
93 // Output slot #1 writes into a TH1 container
94 //DefineOutput(0,TTree::Class());
95 DefineOutput(1,TList::Class());
96
97 // user can change as needed
98 debug = 0;
99 doMC = false;
100 trigger = AliVEvent::kMB;
101 dopupcut = true;
102 doeffcorr = true;
103 centhandle = 1;
104 fbit = 128;
105 zvtxcut = 10.0;
106 centcut = 5.0;
107 nclscut = 60;
108 dcacutz = 5.0;
109 dcacutxy = 5.0;
110 dodcacuts = false;
111 outeta = 0.8;
112 ineta = 0.5;
113 excleta = 0.0;
114 ptmin = 0.2;
115 ptmax = 5.0;
116 doacuts = false;
117 nspid = 2.0;
118 cbinlo = 3;
119 cbinhi = 4;
120 donested = true;
121 dopaircut = true;
122 centlo = 20.0;
123 centhi = 60.0;
124
125}
126
127
128
129//---------------------------------------------------
130AliAnalysisTaskCMEv2A::~AliAnalysisTaskCMEv2A()
131{
132 // Default class destructor
133
134 cout<<"Default class destructor called. Analysis fun time has ended"<<endl;
135
136}
137
138
139
140//---------------------------------------------------
141void AliAnalysisTaskCMEv2A::UserCreateOutputObjects()
142{
143
144 // Create output objects, like histograms - called once
145
146 cout<<"UserCreateOutputObjects called, now making histrograms and things"<<endl;
147
148
149
150 // -------------------------- //
151 // --- create output list --- //
152 // -------------------------- //
153
154 fOutputList = new TList();
155 fOutputList->SetName(GetName());
156 fOutputList->SetOwner(kTRUE);
157
158
159
160 // ------------------------- //
161 // --- create histograms --- //
162 // ------------------------- //
163
164 fHistPt = new TH1F("fHistPt","",50,0.0,5.0);
165 fHistPhi = new TH1F("fHistPhi","",63,0.0,6.3);
166 fHistEta = new TH1F("fHistEta","",300,-1.5,1.5);
167 fHistCharge = new TH1F("fHistCharge","",3,-1.5,1.5);
168 fHistTPCncls = new TH1F("fHistTPCncls","",161,-0.5,160.5);
169 fHistDedx = new TH1F("fHistDedx","",1000,-20,980);
170 fHistDCAxy = new TH1F("fHistDCAxy","",350,-3.5,3.5);
171 fHistDCAz = new TH1F("fHistDCAz","",350,-3.5,3.5);
172 fHistDCAxyAfter = new TH1F("fHistDCAxyAfter","",350,-3.5,3.5);
173 fHistDCAzAfter = new TH1F("fHistDCAzAfter","",350,-3.5,3.5);
174 fOutputList->Add(fHistPt);
175 fOutputList->Add(fHistPhi);
176 fOutputList->Add(fHistEta);
177 fOutputList->Add(fHistCharge);
178 fOutputList->Add(fHistTPCncls);
179 fOutputList->Add(fHistDedx);
180 fOutputList->Add(fHistDCAxy);
181 fOutputList->Add(fHistDCAz);
182 fOutputList->Add(fHistDCAxyAfter);
183 fOutputList->Add(fHistDCAzAfter);
184
185 fHistPosPt = new TH1F("fHistPosPt","",50,0.0,5.0);
186 fHistPosPhi = new TH1F("fHistPosPhi","",63,0.0,6.3);
187 fHistPosEta = new TH1F("fHistPosEta","",300,-1.5,1.5);
188 fOutputList->Add(fHistPosPt);
189 fOutputList->Add(fHistPosPhi);
190 fOutputList->Add(fHistPosEta);
191
192 fHistNegPt = new TH1F("fHistNegPt","",50,0.0,5.0);
193 fHistNegPhi = new TH1F("fHistNegPhi","",63,0.0,6.3);
194 fHistNegEta = new TH1F("fHistNegEta","",300,-1.5,1.5);
195 fOutputList->Add(fHistNegPt);
196 fOutputList->Add(fHistNegPhi);
197 fOutputList->Add(fHistNegEta);
198
199
200
201 fHistPtPion = new TH1F("fHistPtPion","",50,0.0,5.0);
202 fHistPtPionH = new TH1F("fHistPtPionH","",50,0.0,5.0);
203 fHistPtPionL = new TH1F("fHistPtPionL","",50,0.0,5.0);
204 fHistNsigmaPion = new TH1F("fHistNsigmaPion","",100,-5.0,5.0);
205 fHistPtProt = new TH1F("fHistPtProt","",50,0.0,5.0);
206 fHistPtProtH = new TH1F("fHistPtProtH","",50,0.0,5.0);
207 fHistPtProtL = new TH1F("fHistPtProtL","",50,0.0,5.0);
208 fHistNsigmaProt = new TH1F("fHistNsigmaProt","",100,-5.0,5.0);
209 fOutputList->Add(fHistPtPion);
210 fOutputList->Add(fHistPtPionH);
211 fOutputList->Add(fHistPtPionL);
212 fOutputList->Add(fHistNsigmaPion);
213 fOutputList->Add(fHistPtProt);
214 fOutputList->Add(fHistPtProtH);
215 fOutputList->Add(fHistPtProtL);
216 fOutputList->Add(fHistNsigmaProt);
217
218
219 fHistPtMC = new TH1F("fHistPtMC","",50,0.0,5.0);
220 fHistPhiMC = new TH1F("fHistPhiMC","",63,0.0,6.3);
221 fHistEtaMC = new TH1F("fHistEtaMC","",110,-1.1,1.1);
222 fHistChargeMC = new TH1F("fHistChargeMC","",3,-1.5,1.5);
223 fHistTPCnclsMC = new TH1F("fHistTPCnclsMC","",161,-0.5,160.5);
224 fHistDedxMC = new TH1F("fHistDedxMC","",1000,-20,980);
225 fHistDCAxyMC = new TH1F("fHistDCAxyMC","",350,-3.5,3.5);
226 fHistDCAzMC = new TH1F("fHistDCAzMC","",350,-3.5,3.5);
227 fOutputList->Add(fHistPtMC);
228 fOutputList->Add(fHistPhiMC);
229 fOutputList->Add(fHistEtaMC);
230 fOutputList->Add(fHistChargeMC);
231 fOutputList->Add(fHistTPCnclsMC);
232 fOutputList->Add(fHistDedxMC);
233 fOutputList->Add(fHistDCAxyMC);
234 fOutputList->Add(fHistDCAzMC);
235
236
237
238
239 fHistPlaneV0h2 = new TH1F("fHistPlaneV0h2","",640,-3.2,3.2);
240 fHistPlaneV0Ah2 = new TH1F("fHistPlaneV0Ah2","",640,-3.2,3.2);
241 fHistPlaneV0Ch2 = new TH1F("fHistPlaneV0Ch2","",640,-3.2,3.2);
242 fHistPlaneV0ACDCh2 = new TH1F("fHistPlaneV0ACDCh2","",640,-3.2,3.2);
243 fOutputList->Add(fHistPlaneV0h2);
244 fOutputList->Add(fHistPlaneV0Ah2);
245 fOutputList->Add(fHistPlaneV0Ch2);
246 fOutputList->Add(fHistPlaneV0ACDCh2);
247
248
249 fHistZVtx = new TH1F("fHistZVtx","",80,-20,20);
250 fHistZVtxD = new TH1F("fHistZVtxD","",100,-5,5);
251 fHistCentTRK = new TH1F("fHistCentTRK","",100,0,100);
252 fHistCentV0M = new TH1F("fHistCentV0M","",100,0,100);
253 fHistCentDIFF = new TH1F("fHistCentDIFF","",162,-40.5,40.5);
254 fHistCentDIAG = new TH1F("fHistCentDIAG","",6,-3.5,2.5);
255 fOutputList->Add(fHistZVtx);
256 fOutputList->Add(fHistZVtxD);
257 fOutputList->Add(fHistCentTRK);
258 fOutputList->Add(fHistCentV0M);
259 fOutputList->Add(fHistCentDIFF);
260 fOutputList->Add(fHistCentDIAG);
261
262 fHistZVtxMC = new TH1F("fHistZVtxMC","",80,-20,20);
263 fHistZVtxDiffMC = new TH1F("fHistZVtxDiffMC","",100,-5,5);
264 fOutputList->Add(fHistZVtxMC);
265 fOutputList->Add(fHistZVtxDiffMC);
266
267 fHistCtrkDIAG = new TH1F("fHistCtrkDIAG","",1001,-0.5,1000.5);
268 fHistVtxRDIAG = new TH1F("fHistVtxRDIAG","",100,-5,5);
269 fHistVtxSDIAG = new TH1F("fHistVtxSDIAG","",100,-5,5);
270 fHistVtxRDIBG = new TH1F("fHistVtxRDIBG","",100,-5,5);
271 fHistVtxSDIBG = new TH1F("fHistVtxSDIBG","",100,-5,5);
272 fOutputList->Add(fHistCtrkDIAG);
273 fOutputList->Add(fHistVtxRDIAG);
274 fOutputList->Add(fHistVtxSDIAG);
275 fOutputList->Add(fHistVtxRDIBG);
276 fOutputList->Add(fHistVtxSDIBG);
277
278
279 fHistCentTRKAVEkMB = new TH1F("fHistCentTRKAVEkMB","",100,0,100);
280 fHistCentV0MAVEkMB = new TH1F("fHistCentV0MAVEkMB","",100,0,100);
281 fHistCentTRKAVEkCentral = new TH1F("fHistCentTRKAVEkCentral","",100,0,100);
282 fHistCentV0MAVEkCentral = new TH1F("fHistCentV0MAVEkCentral","",100,0,100);
283 fHistCentTRKAVEkSemiCentral = new TH1F("fHistCentTRKAVEkSemiCentral","",100,0,100);
284 fHistCentV0MAVEkSemiCentral = new TH1F("fHistCentV0MAVEkSemiCentral","",100,0,100);
285 fHistCentTRKAVEkA3 = new TH1F("fHistCentTRKAVEkA3","",100,0,100);
286 fHistCentV0MAVEkA3 = new TH1F("fHistCentV0MAVEkA3","",100,0,100);
287 fHistCentTRKAVEkSel = new TH1F("fHistCentTRKAVEkSel","",100,0,100);
288 fHistCentV0MAVEkSel = new TH1F("fHistCentV0MAVEkSel","",100,0,100);
289 fOutputList->Add(fHistCentTRKAVEkMB);
290 fOutputList->Add(fHistCentV0MAVEkMB);
291 fOutputList->Add(fHistCentTRKAVEkCentral);
292 fOutputList->Add(fHistCentV0MAVEkCentral);
293 fOutputList->Add(fHistCentTRKAVEkSemiCentral);
294 fOutputList->Add(fHistCentV0MAVEkSemiCentral);
295 fOutputList->Add(fHistCentTRKAVEkA3);
296 fOutputList->Add(fHistCentV0MAVEkA3);
297 fOutputList->Add(fHistCentTRKAVEkSel);
298 fOutputList->Add(fHistCentV0MAVEkSel);
299
300
301 // -------------------------- //
302 // --- flow and asymmetry --- //
303 // -------------------------- //
304
305 float tpmax = 1e10;
306 float tpmin = -1e10;
307
308 fProfMeanChargePt = new TProfile("fProfMeanChargePt","",50,0,5.0,tpmin,tpmax,"");
309 fProfMeanChargeEta = new TProfile("fProfMeanChargeEta","",300,-1.5,1.5,tpmin,tpmax,"");
310
311 h_eta_pos_F1 = new TH1F("h_eta_pos_F1","",300,-1.5,1.5);
312 h_eta_neg_F1 = new TH1F("h_eta_neg_F1","",300,-1.5,1.5);
313 h_eta_pos_F3 = new TH1F("h_eta_pos_F3","",300,-1.5,1.5);
314 h_eta_neg_F3 = new TH1F("h_eta_neg_F3","",300,-1.5,1.5);
315 fOutputList->Add(h_eta_pos_F1);
316 fOutputList->Add(h_eta_neg_F1);
317 fOutputList->Add(h_eta_pos_F3);
318 fOutputList->Add(h_eta_neg_F3);
319 h_cut_eta_pos_F1 = new TH1F("h_cut_eta_pos_F1","",300,-1.5,1.5);
320 h_cut_eta_neg_F1 = new TH1F("h_cut_eta_neg_F1","",300,-1.5,1.5);
321 h_cut_eta_pos_F3 = new TH1F("h_cut_eta_pos_F3","",300,-1.5,1.5);
322 h_cut_eta_neg_F3 = new TH1F("h_cut_eta_neg_F3","",300,-1.5,1.5);
323 fOutputList->Add(h_cut_eta_pos_F1);
324 fOutputList->Add(h_cut_eta_neg_F1);
325 fOutputList->Add(h_cut_eta_pos_F3);
326 fOutputList->Add(h_cut_eta_neg_F3);
327
328 // --- charge asymmetry (average charge in event)
329 h_AT_eta = new TProfile("h_AT_eta","",150,-1.5,1.5,tpmin,tpmax,"");
330 h_AT_etaF1 = new TProfile("h_AT_etaF1","",150,-1.5,1.5,tpmin,tpmax,"");
331 h_AT_etaF3 = new TProfile("h_AT_etaF3","",150,-1.5,1.5,tpmin,tpmax,"");
332 h_AT_etaMC = new TProfile("h_AT_etaMC","",150,-1.5,1.5,tpmin,tpmax,"");
333 fOutputList->Add(h_AT_eta);
334 fOutputList->Add(h_AT_etaF1);
335 fOutputList->Add(h_AT_etaF3);
336 fOutputList->Add(h_AT_etaMC);
337
338 // --- charge asymmetry (average charge in event)
339 h2_AT_eta = new TH2F("h2_AT_eta","",150,-1.5,1.5,100,-1.0,1.0);
340 h2_AT_etaF1 = new TH2F("h2_AT_etaF1","",150,-1.5,1.5,100,-1.0,1.0);
341 h2_AT_etaF3 = new TH2F("h2_AT_etaF3","",150,-1.5,1.5,100,-1.0,1.0);
342 h2_AT_etaMC = new TH2F("h2_AT_etaMC","",150,-1.5,1.5,100,-1.0,1.0);
343 fOutputList->Add(h2_AT_eta);
344 fOutputList->Add(h2_AT_etaF1);
345 fOutputList->Add(h2_AT_etaF3);
346 fOutputList->Add(h2_AT_etaMC);
347
348 // --- charge asymmetry (average charge in event)
349 h_AT_cut_eta = new TProfile("h_AT_cut_eta","",150,-1.5,1.5,tpmin,tpmax,"");
350 h_AT_cut_etaF1 = new TProfile("h_AT_cut_etaF1","",150,-1.5,1.5,tpmin,tpmax,"");
351 h_AT_cut_etaF3 = new TProfile("h_AT_cut_etaF3","",150,-1.5,1.5,tpmin,tpmax,"");
352 h_AT_cut_etaMC = new TProfile("h_AT_cut_etaMC","",150,-1.5,1.5,tpmin,tpmax,"");
353 fOutputList->Add(h_AT_cut_eta);
354 fOutputList->Add(h_AT_cut_etaF1);
355 fOutputList->Add(h_AT_cut_etaF3);
356 fOutputList->Add(h_AT_cut_etaMC);
357
358 // --- charge asymmetry (average charge in event)
359 h2_AT_cut_eta = new TH2F("h2_AT_cut_eta","",150,-1.5,1.5,100,-1.0,1.0);
360 h2_AT_cut_etaF1 = new TH2F("h2_AT_cut_etaF1","",150,-1.5,1.5,100,-1.0,1.0);
361 h2_AT_cut_etaF3 = new TH2F("h2_AT_cut_etaF3","",150,-1.5,1.5,100,-1.0,1.0);
362 h2_AT_cut_etaMC = new TH2F("h2_AT_cut_etaMC","",150,-1.5,1.5,100,-1.0,1.0);
363 fOutputList->Add(h2_AT_cut_eta);
364 fOutputList->Add(h2_AT_cut_etaF1);
365 fOutputList->Add(h2_AT_cut_etaF3);
366 fOutputList->Add(h2_AT_cut_etaMC);
367
368 // --- charge asymmetry (average charge in event)
369 h_A_cent = new TProfile("h_A_cent","",100,0,100,tpmin,tpmax,"");
370 h_rA_cent = new TProfile("h_rA_cent","",100,0,100,tpmin,tpmax,"");
371 h_r1A_cent = new TProfile("h_r1A_cent","",100,0,100,tpmin,tpmax,"");
372 h_r2A_cent = new TProfile("h_r2A_cent","",100,0,100,tpmin,tpmax,"");
373 h_A_centF1 = new TProfile("h_A_centF1","",100,0,100,tpmin,tpmax,"");
374 h_A_centF3 = new TProfile("h_A_centF3","",100,0,100,tpmin,tpmax,"");
375 h_A_centMC = new TProfile("h_A_centMC","",100,0,100,tpmin,tpmax,"");
376 fOutputList->Add(h_A_cent);
377 fOutputList->Add(h_rA_cent);
378 fOutputList->Add(h_r1A_cent);
379 fOutputList->Add(h_r2A_cent);
380 fOutputList->Add(h_A_centF1);
381 fOutputList->Add(h_A_centF3);
382 fOutputList->Add(h_A_centMC);
383
384 // --- charge asymmetry (average charge in event)
385 h2_A_cent = new TH2F("h2_A_cent","",100,0,100,100,-1.0,1.0);
386 h2_rA_cent = new TH2F("h2_rA_cent","",100,0,100,100,-1.0,1.0);
387 h2_r1A_cent = new TH2F("h2_r1A_cent","",100,0,100,100,-1.0,1.0);
388 h2_r2A_cent = new TH2F("h2_r2A_cent","",100,0,100,100,-1.0,1.0);
389 h2_A_centF1 = new TH2F("h2_A_centF1","",100,0,100,100,-1.0,1.0);
390 h2_A_centF3 = new TH2F("h2_A_centF3","",100,0,100,100,-1.0,1.0);
391 h2_A_centMC = new TH2F("h2_A_centMC","",100,0,100,100,-1.0,1.0);
392 fOutputList->Add(h2_A_cent);
393 fOutputList->Add(h2_rA_cent);
394 fOutputList->Add(h2_r1A_cent);
395 fOutputList->Add(h2_r2A_cent);
396 fOutputList->Add(h2_A_centF1);
397 fOutputList->Add(h2_A_centF3);
398 fOutputList->Add(h2_A_centMC);
399
400 // --- MONTE CARLO second harmonic vs centrality
401 h_MCq22_cent = new TProfile("h_MCq22_cent","",100,0,100,tpmin,tpmax,"");
402 h_MCq22_centP = new TProfile("h_MCq22_centP","",100,0,100,tpmin,tpmax,"");
403 h_MCq22_centN = new TProfile("h_MCq22_centN","",100,0,100,tpmin,tpmax,"");
404 fOutputList->Add(h_MCq22_cent);
405 fOutputList->Add(h_MCq22_centP);
406 fOutputList->Add(h_MCq22_centN);
407 h_MCAq22_cent = new TProfile("h_MCAq22_cent","",100,0,100,tpmin,tpmax,"");
408 h_MCAq22_centP = new TProfile("h_MCAq22_centP","",100,0,100,tpmin,tpmax,"");
409 h_MCAq22_centN = new TProfile("h_MCAq22_centN","",100,0,100,tpmin,tpmax,"");
410 fOutputList->Add(h_MCAq22_cent);
411 fOutputList->Add(h_MCAq22_centP);
412 fOutputList->Add(h_MCAq22_centN);
413 h_MCq22gap0_cent = new TProfile("h_MCq22gap0_cent","",100,0,100,tpmin,tpmax,"");
414 h_MCq22gap0_centP = new TProfile("h_MCq22gap0_centP","",100,0,100,tpmin,tpmax,"");
415 h_MCq22gap0_centN = new TProfile("h_MCq22gap0_centN","",100,0,100,tpmin,tpmax,"");
416 fOutputList->Add(h_MCq22gap0_cent);
417 fOutputList->Add(h_MCq22gap0_centP);
418 fOutputList->Add(h_MCq22gap0_centN);
419 h_MCq22gap1_cent = new TProfile("h_MCq22gap1_cent","",100,0,100,tpmin,tpmax,"");
420 h_MCq22gap1_centP = new TProfile("h_MCq22gap1_centP","",100,0,100,tpmin,tpmax,"");
421 h_MCq22gap1_centN = new TProfile("h_MCq22gap1_centN","",100,0,100,tpmin,tpmax,"");
422 fOutputList->Add(h_MCq22gap1_cent);
423 fOutputList->Add(h_MCq22gap1_centP);
424 fOutputList->Add(h_MCq22gap1_centN);
425
426 hMC_AT_X_deta = new TProfile("hMC_AT_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
427 hMC_AT_X_detaP = new TProfile("hMC_AT_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
428 hMC_AT_X_detaN = new TProfile("hMC_AT_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
429 fOutputList->Add(hMC_AT_X_deta);
430 fOutputList->Add(hMC_AT_X_detaP);
431 fOutputList->Add(hMC_AT_X_detaN);
432 //hMC_diffq22_X_deta = new TProfile("hMC_diffq22_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
433 hMC_diffq22_X_detaP = new TProfile("hMC_diffq22_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
434 hMC_diffq22_X_detaN = new TProfile("hMC_diffq22_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
435 //hMC_ATdiffq22_X_deta = new TProfile("hMC_ATdiffq22_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
436 hMC_ATdiffq22_X_detaP = new TProfile("hMC_ATdiffq22_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
437 hMC_ATdiffq22_X_detaN = new TProfile("hMC_ATdiffq22_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
438 fOutputList->Add(hMC_diffq22_X_detaP);
439 fOutputList->Add(hMC_diffq22_X_detaN);
440 fOutputList->Add(hMC_ATdiffq22_X_detaP);
441 fOutputList->Add(hMC_ATdiffq22_X_detaN);
442 //hMC_diffq32_X_deta = new TProfile("hMC_diffq32_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
443 hMC_diffq32_X_detaP = new TProfile("hMC_diffq32_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
444 hMC_diffq32_X_detaN = new TProfile("hMC_diffq32_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
445 //hMC_ATdiffq32_X_deta = new TProfile("hMC_ATdiffq32_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
446 hMC_ATdiffq32_X_detaP = new TProfile("hMC_ATdiffq32_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
447 hMC_ATdiffq32_X_detaN = new TProfile("hMC_ATdiffq32_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
448 fOutputList->Add(hMC_diffq32_X_detaP);
449 fOutputList->Add(hMC_diffq32_X_detaN);
450 fOutputList->Add(hMC_ATdiffq32_X_detaP);
451 fOutputList->Add(hMC_ATdiffq32_X_detaN);
452 //hMC_diffq42_X_deta = new TProfile("hMC_diffq42_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
453 hMC_diffq42_X_detaP = new TProfile("hMC_diffq42_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
454 hMC_diffq42_X_detaN = new TProfile("hMC_diffq42_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
455 //hMC_ATdiffq42_X_deta = new TProfile("hMC_ATdiffq42_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
456 hMC_ATdiffq42_X_detaP = new TProfile("hMC_ATdiffq42_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
457 hMC_ATdiffq42_X_detaN = new TProfile("hMC_ATdiffq42_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
458 fOutputList->Add(hMC_diffq42_X_detaP);
459 fOutputList->Add(hMC_diffq42_X_detaN);
460 fOutputList->Add(hMC_ATdiffq42_X_detaP);
461 fOutputList->Add(hMC_ATdiffq42_X_detaN);
462
463
464 // --- q-vector components to check recentering
465 h_X22_cent = new TProfile("h_X22_cent","",100,0,100,tpmin,tpmax,"");
466 h_X22_centP = new TProfile("h_X22_centP","",100,0,100,tpmin,tpmax,"");
467 h_X22_centN = new TProfile("h_X22_centN","",100,0,100,tpmin,tpmax,"");
468 h_Y22_cent = new TProfile("h_Y22_cent","",100,0,100,tpmin,tpmax,"");
469 h_Y22_centP = new TProfile("h_Y22_centP","",100,0,100,tpmin,tpmax,"");
470 h_Y22_centN = new TProfile("h_Y22_centN","",100,0,100,tpmin,tpmax,"");
471 fOutputList->Add(h_X22_cent);
472 fOutputList->Add(h_X22_centP);
473 fOutputList->Add(h_X22_centN);
474 fOutputList->Add(h_Y22_cent);
475 fOutputList->Add(h_Y22_centP);
476 fOutputList->Add(h_Y22_centN);
477 // --- outer left
478 h_X22_lo_cent = new TProfile("h_X22_lo_cent","",100,0,100,tpmin,tpmax,"");
479 h_X22_lo_centP = new TProfile("h_X22_lo_centP","",100,0,100,tpmin,tpmax,"");
480 h_X22_lo_centN = new TProfile("h_X22_lo_centN","",100,0,100,tpmin,tpmax,"");
481 h_Y22_lo_cent = new TProfile("h_Y22_lo_cent","",100,0,100,tpmin,tpmax,"");
482 h_Y22_lo_centP = new TProfile("h_Y22_lo_centP","",100,0,100,tpmin,tpmax,"");
483 h_Y22_lo_centN = new TProfile("h_Y22_lo_centN","",100,0,100,tpmin,tpmax,"");
484 fOutputList->Add(h_X22_lo_cent);
485 fOutputList->Add(h_X22_lo_centP);
486 fOutputList->Add(h_X22_lo_centN);
487 fOutputList->Add(h_Y22_lo_cent);
488 fOutputList->Add(h_Y22_lo_centP);
489 fOutputList->Add(h_Y22_lo_centN);
490 // --- inner left
491 h_X22_li_cent = new TProfile("h_X22_li_cent","",100,0,100,tpmin,tpmax,"");
492 h_X22_li_centP = new TProfile("h_X22_li_centP","",100,0,100,tpmin,tpmax,"");
493 h_X22_li_centN = new TProfile("h_X22_li_centN","",100,0,100,tpmin,tpmax,"");
494 h_Y22_li_cent = new TProfile("h_Y22_li_cent","",100,0,100,tpmin,tpmax,"");
495 h_Y22_li_centP = new TProfile("h_Y22_li_centP","",100,0,100,tpmin,tpmax,"");
496 h_Y22_li_centN = new TProfile("h_Y22_li_centN","",100,0,100,tpmin,tpmax,"");
497 fOutputList->Add(h_X22_li_cent);
498 fOutputList->Add(h_X22_li_centP);
499 fOutputList->Add(h_X22_li_centN);
500 fOutputList->Add(h_Y22_li_cent);
501 fOutputList->Add(h_Y22_li_centP);
502 fOutputList->Add(h_Y22_li_centN);
503 // --- inner right
504 h_X22_ri_cent = new TProfile("h_X22_ri_cent","",100,0,100,tpmin,tpmax,"");
505 h_X22_ri_centP = new TProfile("h_X22_ri_centP","",100,0,100,tpmin,tpmax,"");
506 h_X22_ri_centN = new TProfile("h_X22_ri_centN","",100,0,100,tpmin,tpmax,"");
507 h_Y22_ri_cent = new TProfile("h_Y22_ri_cent","",100,0,100,tpmin,tpmax,"");
508 h_Y22_ri_centP = new TProfile("h_Y22_ri_centP","",100,0,100,tpmin,tpmax,"");
509 h_Y22_ri_centN = new TProfile("h_Y22_ri_centN","",100,0,100,tpmin,tpmax,"");
510 fOutputList->Add(h_X22_ri_cent);
511 fOutputList->Add(h_X22_ri_centP);
512 fOutputList->Add(h_X22_ri_centN);
513 fOutputList->Add(h_Y22_ri_cent);
514 fOutputList->Add(h_Y22_ri_centP);
515 fOutputList->Add(h_Y22_ri_centN);
516 // --- outer right
517 h_X22_ro_cent = new TProfile("h_X22_ro_cent","",100,0,100,tpmin,tpmax,"");
518 h_X22_ro_centP = new TProfile("h_X22_ro_centP","",100,0,100,tpmin,tpmax,"");
519 h_X22_ro_centN = new TProfile("h_X22_ro_centN","",100,0,100,tpmin,tpmax,"");
520 h_Y22_ro_cent = new TProfile("h_Y22_ro_cent","",100,0,100,tpmin,tpmax,"");
521 h_Y22_ro_centP = new TProfile("h_Y22_ro_centP","",100,0,100,tpmin,tpmax,"");
522 h_Y22_ro_centN = new TProfile("h_Y22_ro_centN","",100,0,100,tpmin,tpmax,"");
523 fOutputList->Add(h_X22_ro_cent);
524 fOutputList->Add(h_X22_ro_centP);
525 fOutputList->Add(h_X22_ro_centN);
526 fOutputList->Add(h_Y22_ro_cent);
527 fOutputList->Add(h_Y22_ro_centP);
528 fOutputList->Add(h_Y22_ro_centN);
529
530 // --- second harmonic vs centrality
531 h_q22_cent = new TProfile("h_q22_cent","",100,0,100,tpmin,tpmax,"");
532 h_q22_centP = new TProfile("h_q22_centP","",100,0,100,tpmin,tpmax,"");
533 h_q22_centN = new TProfile("h_q22_centN","",100,0,100,tpmin,tpmax,"");
534 fOutputList->Add(h_q22_cent);
535 fOutputList->Add(h_q22_centP);
536 fOutputList->Add(h_q22_centN);
537 h_q23_cent = new TProfile("h_q23_cent","",100,0,100,tpmin,tpmax,"");
538 h_q23_centP = new TProfile("h_q23_centP","",100,0,100,tpmin,tpmax,"");
539 h_q23_centN = new TProfile("h_q23_centN","",100,0,100,tpmin,tpmax,"");
540 fOutputList->Add(h_q23_cent);
541 fOutputList->Add(h_q23_centP);
542 fOutputList->Add(h_q23_centN);
543 h_q24_cent = new TProfile("h_q24_cent","",100,0,100,tpmin,tpmax,"");
544 h_q24_centP = new TProfile("h_q24_centP","",100,0,100,tpmin,tpmax,"");
545 h_q24_centN = new TProfile("h_q24_centN","",100,0,100,tpmin,tpmax,"");
546 fOutputList->Add(h_q24_cent);
547 fOutputList->Add(h_q24_centP);
548 fOutputList->Add(h_q24_centN);
549 h_q22gap0_cent = new TProfile("h_q22gap0_cent","",100,0,100,tpmin,tpmax,"");
550 h_q22gap0_centP = new TProfile("h_q22gap0_centP","",100,0,100,tpmin,tpmax,"");
551 h_q22gap0_centN = new TProfile("h_q22gap0_centN","",100,0,100,tpmin,tpmax,"");
552 fOutputList->Add(h_q22gap0_cent);
553 fOutputList->Add(h_q22gap0_centP);
554 fOutputList->Add(h_q22gap0_centN);
555 h_q22gap1_cent = new TProfile("h_q22gap1_cent","",100,0,100,tpmin,tpmax,"");
556 h_q22gap1_centP = new TProfile("h_q22gap1_centP","",100,0,100,tpmin,tpmax,"");
557 h_q22gap1_centN = new TProfile("h_q22gap1_centN","",100,0,100,tpmin,tpmax,"");
558 fOutputList->Add(h_q22gap1_cent);
559 fOutputList->Add(h_q22gap1_centP);
560 fOutputList->Add(h_q22gap1_centN);
561
562 // --- third harmonic vs centrality
563 h_q32_cent = new TProfile("h_q32_cent","",100,0,100,tpmin,tpmax,"");
564 h_q32_centP = new TProfile("h_q32_centP","",100,0,100,tpmin,tpmax,"");
565 h_q32_centN = new TProfile("h_q32_centN","",100,0,100,tpmin,tpmax,"");
566 fOutputList->Add(h_q32_cent);
567 fOutputList->Add(h_q32_centP);
568 fOutputList->Add(h_q32_centN);
569 h_q32gap0_cent = new TProfile("h_q32gap0_cent","",100,0,100,tpmin,tpmax,"");
570 h_q32gap0_centP = new TProfile("h_q32gap0_centP","",100,0,100,tpmin,tpmax,"");
571 h_q32gap0_centN = new TProfile("h_q32gap0_centN","",100,0,100,tpmin,tpmax,"");
572 fOutputList->Add(h_q32gap0_cent);
573 fOutputList->Add(h_q32gap0_centP);
574 fOutputList->Add(h_q32gap0_centN);
575 h_q32gap1_cent = new TProfile("h_q32gap1_cent","",100,0,100,tpmin,tpmax,"");
576 h_q32gap1_centP = new TProfile("h_q32gap1_centP","",100,0,100,tpmin,tpmax,"");
577 h_q32gap1_centN = new TProfile("h_q32gap1_centN","",100,0,100,tpmin,tpmax,"");
578 fOutputList->Add(h_q32gap1_cent);
579 fOutputList->Add(h_q32gap1_centP);
580 fOutputList->Add(h_q32gap1_centN);
581
582 // --- fourth harmonic vs centrality
583 h_q42_cent = new TProfile("h_q42_cent","",100,0,100,tpmin,tpmax,"");
584 h_q42_centP = new TProfile("h_q42_centP","",100,0,100,tpmin,tpmax,"");
585 h_q42_centN = new TProfile("h_q42_centN","",100,0,100,tpmin,tpmax,"");
586 fOutputList->Add(h_q42_cent);
587 fOutputList->Add(h_q42_centP);
588 fOutputList->Add(h_q42_centN);
589 h_q42gap0_cent = new TProfile("h_q42gap0_cent","",100,0,100,tpmin,tpmax,"");
590 h_q42gap0_centP = new TProfile("h_q42gap0_centP","",100,0,100,tpmin,tpmax,"");
591 h_q42gap0_centN = new TProfile("h_q42gap0_centN","",100,0,100,tpmin,tpmax,"");
592 fOutputList->Add(h_q42gap0_cent);
593 fOutputList->Add(h_q42gap0_centP);
594 fOutputList->Add(h_q42gap0_centN);
595 h_q42gap1_cent = new TProfile("h_q42gap1_cent","",100,0,100,tpmin,tpmax,"");
596 h_q42gap1_centP = new TProfile("h_q42gap1_centP","",100,0,100,tpmin,tpmax,"");
597 h_q42gap1_centN = new TProfile("h_q42gap1_centN","",100,0,100,tpmin,tpmax,"");
598 fOutputList->Add(h_q42gap1_cent);
599 fOutputList->Add(h_q42gap1_centP);
600 fOutputList->Add(h_q42gap1_centN);
601
602 // harmonics vs charge asymmetry in centrality bins
603 for(int i=0; i<10; i++)
604 {
605 // --- second harmonic
606 h_q22qasym_cent[i] = new TProfile(Form("h_q22qasym_cent_%d",i),Form("h_q22qasym_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
607 h_q22qasymP_cent[i] = new TProfile(Form("h_q22qasymP_cent_%d",i),Form("h_q22qasymP_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
608 h_q22qasymN_cent[i] = new TProfile(Form("h_q22qasymN_cent_%d",i),Form("h_q22qasymN_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
609 h_q22qasym_gap0_cent[i] = new TProfile(Form("h_q22qasym_gap0_cent_%d",i),Form("h_q22qasym_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
610 h_q22qasymP_gap0_cent[i] = new TProfile(Form("h_q22qasymP_gap0_cent_%d",i),Form("h_q22qasymP_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
611 h_q22qasymN_gap0_cent[i] = new TProfile(Form("h_q22qasymN_gap0_cent_%d",i),Form("h_q22qasymN_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
612 h_q22qasym_gap1_cent[i] = new TProfile(Form("h_q22qasym_gap1_cent_%d",i),Form("h_q22qasym_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
613 h_q22qasymP_gap1_cent[i] = new TProfile(Form("h_q22qasymP_gap1_cent_%d",i),Form("h_q22qasymP_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
614 h_q22qasymN_gap1_cent[i] = new TProfile(Form("h_q22qasymN_gap1_cent_%d",i),Form("h_q22qasymN_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
615 h_q24qasym_cent[i] = new TProfile(Form("h_q24qasym_cent_%d",i),Form("h_q24qasym_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
616 h_q24qasymP_cent[i] = new TProfile(Form("h_q24qasymP_cent_%d",i),Form("h_q24qasymP_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
617 h_q24qasymN_cent[i] = new TProfile(Form("h_q24qasymN_cent_%d",i),Form("h_q24qasymN_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
618
619 // --- third harmonic
620 h_q32qasym_cent[i] = new TProfile(Form("h_q32qasym_cent_%d",i),Form("h_q32qasym_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
621 h_q32qasymP_cent[i] = new TProfile(Form("h_q32qasymP_cent_%d",i),Form("h_q32qasymP_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
622 h_q32qasymN_cent[i] = new TProfile(Form("h_q32qasymN_cent_%d",i),Form("h_q32qasymN_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
623 h_q32qasym_gap0_cent[i] = new TProfile(Form("h_q32qasym_gap0_cent_%d",i),Form("h_q32qasym_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
624 h_q32qasymP_gap0_cent[i] = new TProfile(Form("h_q32qasymP_gap0_cent_%d",i),Form("h_q32qasymP_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
625 h_q32qasymN_gap0_cent[i] = new TProfile(Form("h_q32qasymN_gap0_cent_%d",i),Form("h_q32qasymN_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
626 h_q32qasym_gap1_cent[i] = new TProfile(Form("h_q32qasym_gap1_cent_%d",i),Form("h_q32qasym_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
627 h_q32qasymP_gap1_cent[i] = new TProfile(Form("h_q32qasymP_gap1_cent_%d",i),Form("h_q32qasymP_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
628 h_q32qasymN_gap1_cent[i] = new TProfile(Form("h_q32qasymN_gap1_cent_%d",i),Form("h_q32qasymN_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
629 h_q34qasym_cent[i] = new TProfile(Form("h_q34qasym_cent_%d",i),Form("h_q34qasym_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
630 h_q34qasymP_cent[i] = new TProfile(Form("h_q34qasymP_cent_%d",i),Form("h_q34qasymP_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
631 h_q34qasymN_cent[i] = new TProfile(Form("h_q34qasymN_cent_%d",i),Form("h_q34qasymN_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
632
633 // --- fourth harmonic
634 h_q42qasym_cent[i] = new TProfile(Form("h_q42qasym_cent_%d",i),Form("h_q42qasym_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
635 h_q42qasymP_cent[i] = new TProfile(Form("h_q42qasymP_cent_%d",i),Form("h_q42qasymP_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
636 h_q42qasymN_cent[i] = new TProfile(Form("h_q42qasymN_cent_%d",i),Form("h_q42qasymN_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
637 h_q42qasym_gap0_cent[i] = new TProfile(Form("h_q42qasym_gap0_cent_%d",i),Form("h_q42qasym_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
638 h_q42qasymP_gap0_cent[i] = new TProfile(Form("h_q42qasymP_gap0_cent_%d",i),Form("h_q42qasymP_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
639 h_q42qasymN_gap0_cent[i] = new TProfile(Form("h_q42qasymN_gap0_cent_%d",i),Form("h_q42qasymN_gap0_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
640 h_q42qasym_gap1_cent[i] = new TProfile(Form("h_q42qasym_gap1_cent_%d",i),Form("h_q42qasym_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
641 h_q42qasymP_gap1_cent[i] = new TProfile(Form("h_q42qasymP_gap1_cent_%d",i),Form("h_q42qasymP_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
642 h_q42qasymN_gap1_cent[i] = new TProfile(Form("h_q42qasymN_gap1_cent_%d",i),Form("h_q42qasymN_gap1_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
643 h_q44qasym_cent[i] = new TProfile(Form("h_q44qasym_cent_%d",i),Form("h_q44qasym_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
644 h_q44qasymP_cent[i] = new TProfile(Form("h_q44qasymP_cent_%d",i),Form("h_q44qasymP_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
645 h_q44qasymN_cent[i] = new TProfile(Form("h_q44qasymN_cent_%d",i),Form("h_q44qasymN_cent_%d",i),220,-1.1,1.1,tpmin,tpmax,"");
646 }
647
648
649 // harmonics vs charge asymmetry in centrality bins
650 for(int i=cbinlo; i<cbinhi; i++)
651 {
652 // --- second harmonic
653 fOutputList->Add(h_q22qasym_cent[i]);
654 fOutputList->Add(h_q22qasymP_cent[i]);
655 fOutputList->Add(h_q22qasymN_cent[i]);
656 fOutputList->Add(h_q22qasym_gap0_cent[i]);
657 fOutputList->Add(h_q22qasymP_gap0_cent[i]);
658 fOutputList->Add(h_q22qasymN_gap0_cent[i]);
659 fOutputList->Add(h_q22qasym_gap1_cent[i]);
660 fOutputList->Add(h_q22qasymP_gap1_cent[i]);
661 fOutputList->Add(h_q22qasymN_gap1_cent[i]);
662 fOutputList->Add(h_q24qasym_cent[i]);
663 fOutputList->Add(h_q24qasymP_cent[i]);
664 fOutputList->Add(h_q24qasymN_cent[i]);
665
666 // --- third harmonic
667 fOutputList->Add(h_q32qasym_cent[i]);
668 fOutputList->Add(h_q32qasymP_cent[i]);
669 fOutputList->Add(h_q32qasymN_cent[i]);
670 fOutputList->Add(h_q32qasym_gap0_cent[i]);
671 fOutputList->Add(h_q32qasymP_gap0_cent[i]);
672 fOutputList->Add(h_q32qasymN_gap0_cent[i]);
673 fOutputList->Add(h_q32qasym_gap1_cent[i]);
674 fOutputList->Add(h_q32qasymP_gap1_cent[i]);
675 fOutputList->Add(h_q32qasymN_gap1_cent[i]);
676 fOutputList->Add(h_q34qasym_cent[i]);
677 fOutputList->Add(h_q34qasymP_cent[i]);
678 fOutputList->Add(h_q34qasymN_cent[i]);
679
680 // --- fourth harmonic
681 fOutputList->Add(h_q42qasym_cent[i]);
682 fOutputList->Add(h_q42qasymP_cent[i]);
683 fOutputList->Add(h_q42qasymN_cent[i]);
684 fOutputList->Add(h_q42qasym_gap0_cent[i]);
685 fOutputList->Add(h_q42qasymP_gap0_cent[i]);
686 fOutputList->Add(h_q42qasymN_gap0_cent[i]);
687 fOutputList->Add(h_q42qasym_gap1_cent[i]);
688 fOutputList->Add(h_q42qasymP_gap1_cent[i]);
689 fOutputList->Add(h_q42qasymN_gap1_cent[i]);
690 fOutputList->Add(h_q44qasym_cent[i]);
691 fOutputList->Add(h_q44qasymP_cent[i]);
692 fOutputList->Add(h_q44qasymN_cent[i]);
693 }
694
695
696
697 // --------------------------------------- //
698 // --- three-particle-like correlators --- //
699 // --------------------------------------- //
700
701
702 // --- second harmonic
703 h_Aq22_cent = new TProfile("h_Aq22_cent","",100,0,100,tpmin,tpmax,"");
704 h_Aq22_centP = new TProfile("h_Aq22_centP","",100,0,100,tpmin,tpmax,"");
705 h_Aq22_centN = new TProfile("h_Aq22_centN","",100,0,100,tpmin,tpmax,"");
706 fOutputList->Add(h_Aq22_cent);
707 fOutputList->Add(h_Aq22_centP);
708 fOutputList->Add(h_Aq22_centN);
709 h_Aq22_SL_cent = new TProfile("h_Aq22_SL_cent","",100,0,100,tpmin,tpmax,"");
710 h_Aq22_SL_centP = new TProfile("h_Aq22_SL_centP","",100,0,100,tpmin,tpmax,"");
711 h_Aq22_SL_centN = new TProfile("h_Aq22_SL_centN","",100,0,100,tpmin,tpmax,"");
712 fOutputList->Add(h_Aq22_SL_cent);
713 fOutputList->Add(h_Aq22_SL_centP);
714 fOutputList->Add(h_Aq22_SL_centN);
715 h_Aq22_SLL_cent = new TProfile("h_Aq22_SLL_cent","",100,0,100,tpmin,tpmax,"");
716 h_Aq22_SLL_centP = new TProfile("h_Aq22_SLL_centP","",100,0,100,tpmin,tpmax,"");
717 h_Aq22_SLL_centN = new TProfile("h_Aq22_SLL_centN","",100,0,100,tpmin,tpmax,"");
718 fOutputList->Add(h_Aq22_SLL_cent);
719 fOutputList->Add(h_Aq22_SLL_centP);
720 fOutputList->Add(h_Aq22_SLL_centN);
721 h_Aq22_SLR_cent = new TProfile("h_Aq22_SLR_cent","",100,0,100,tpmin,tpmax,"");
722 h_Aq22_SLR_centP = new TProfile("h_Aq22_SLR_centP","",100,0,100,tpmin,tpmax,"");
723 h_Aq22_SLR_centN = new TProfile("h_Aq22_SLR_centN","",100,0,100,tpmin,tpmax,"");
724 fOutputList->Add(h_Aq22_SLR_cent);
725 fOutputList->Add(h_Aq22_SLR_centP);
726 fOutputList->Add(h_Aq22_SLR_centN);
727
728 // --- third harmonic
729 h_Aq32_cent = new TProfile("h_Aq32_cent","",100,0,100,tpmin,tpmax,"");
730 h_Aq32_centP = new TProfile("h_Aq32_centP","",100,0,100,tpmin,tpmax,"");
731 h_Aq32_centN = new TProfile("h_Aq32_centN","",100,0,100,tpmin,tpmax,"");
732 fOutputList->Add(h_Aq32_cent);
733 fOutputList->Add(h_Aq32_centP);
734 fOutputList->Add(h_Aq32_centN);
735 h_Aq32_SL_cent = new TProfile("h_Aq32_SL_cent","",100,0,100,tpmin,tpmax,"");
736 h_Aq32_SL_centP = new TProfile("h_Aq32_SL_centP","",100,0,100,tpmin,tpmax,"");
737 h_Aq32_SL_centN = new TProfile("h_Aq32_SL_centN","",100,0,100,tpmin,tpmax,"");
738 fOutputList->Add(h_Aq32_SL_cent);
739 fOutputList->Add(h_Aq32_SL_centP);
740 fOutputList->Add(h_Aq32_SL_centN);
741 h_Aq32_SLL_cent = new TProfile("h_Aq32_SLL_cent","",100,0,100,tpmin,tpmax,"");
742 h_Aq32_SLL_centP = new TProfile("h_Aq32_SLL_centP","",100,0,100,tpmin,tpmax,"");
743 h_Aq32_SLL_centN = new TProfile("h_Aq32_SLL_centN","",100,0,100,tpmin,tpmax,"");
744 fOutputList->Add(h_Aq32_SLL_cent);
745 fOutputList->Add(h_Aq32_SLL_centP);
746 fOutputList->Add(h_Aq32_SLL_centN);
747 h_Aq32_SLR_cent = new TProfile("h_Aq32_SLR_cent","",100,0,100,tpmin,tpmax,"");
748 h_Aq32_SLR_centP = new TProfile("h_Aq32_SLR_centP","",100,0,100,tpmin,tpmax,"");
749 h_Aq32_SLR_centN = new TProfile("h_Aq32_SLR_centN","",100,0,100,tpmin,tpmax,"");
750 fOutputList->Add(h_Aq32_SLR_cent);
751 fOutputList->Add(h_Aq32_SLR_centP);
752 fOutputList->Add(h_Aq32_SLR_centN);
753
754 // --- fourth harmonic
755 h_Aq42_cent = new TProfile("h_Aq42_cent","",100,0,100,tpmin,tpmax,"");
756 h_Aq42_centP = new TProfile("h_Aq42_centP","",100,0,100,tpmin,tpmax,"");
757 h_Aq42_centN = new TProfile("h_Aq42_centN","",100,0,100,tpmin,tpmax,"");
758 fOutputList->Add(h_Aq42_cent);
759 fOutputList->Add(h_Aq42_centP);
760 fOutputList->Add(h_Aq42_centN);
761 h_Aq42_SL_cent = new TProfile("h_Aq42_SL_cent","",100,0,100,tpmin,tpmax,"");
762 h_Aq42_SL_centP = new TProfile("h_Aq42_SL_centP","",100,0,100,tpmin,tpmax,"");
763 h_Aq42_SL_centN = new TProfile("h_Aq42_SL_centN","",100,0,100,tpmin,tpmax,"");
764 fOutputList->Add(h_Aq42_SL_cent);
765 fOutputList->Add(h_Aq42_SL_centP);
766 fOutputList->Add(h_Aq42_SL_centN);
767 h_Aq42_SLL_cent = new TProfile("h_Aq42_SLL_cent","",100,0,100,tpmin,tpmax,"");
768 h_Aq42_SLL_centP = new TProfile("h_Aq42_SLL_centP","",100,0,100,tpmin,tpmax,"");
769 h_Aq42_SLL_centN = new TProfile("h_Aq42_SLL_centN","",100,0,100,tpmin,tpmax,"");
770 fOutputList->Add(h_Aq42_SLL_cent);
771 fOutputList->Add(h_Aq42_SLL_centP);
772 fOutputList->Add(h_Aq42_SLL_centN);
773 h_Aq42_SLR_cent = new TProfile("h_Aq42_SLR_cent","",100,0,100,tpmin,tpmax,"");
774 h_Aq42_SLR_centP = new TProfile("h_Aq42_SLR_centP","",100,0,100,tpmin,tpmax,"");
775 h_Aq42_SLR_centN = new TProfile("h_Aq42_SLR_centN","",100,0,100,tpmin,tpmax,"");
776 fOutputList->Add(h_Aq42_SLR_cent);
777 fOutputList->Add(h_Aq42_SLR_centP);
778 fOutputList->Add(h_Aq42_SLR_centN);
779
780
781
782
783 // ------------------------------------
784 // --- differential cumulant histograms
785 // ------------------------------------
786
787 h_diffq22_pt = new TProfile("h_diffq22_pt","",50,0,5,tpmin,tpmax,"");
788 h_diffq22_ptP = new TProfile("h_diffq22_ptP","",50,0,5,tpmin,tpmax,"");
789 h_diffq22_ptN = new TProfile("h_diffq22_ptN","",50,0,5,tpmin,tpmax,"");
790 h_diffAq22_pt = new TProfile("h_diffAq22_pt","",50,0,5,tpmin,tpmax,"");
791 h_diffAq22_ptP = new TProfile("h_diffAq22_ptP","",50,0,5,tpmin,tpmax,"");
792 h_diffAq22_ptN = new TProfile("h_diffAq22_ptN","",50,0,5,tpmin,tpmax,"");
793 fOutputList->Add(h_diffq22_pt);
794 fOutputList->Add(h_diffq22_ptP);
795 fOutputList->Add(h_diffq22_ptN);
796 fOutputList->Add(h_diffAq22_pt);
797 fOutputList->Add(h_diffAq22_ptP);
798 fOutputList->Add(h_diffAq22_ptN);
799 h_diffq22_eta = new TProfile("h_diffq22_eta","",32,-0.8,0.8,tpmin,tpmax,"");
800 h_diffq22_etaP = new TProfile("h_diffq22_etaP","",32,-0.8,0.8,tpmin,tpmax,"");
801 h_diffq22_etaN = new TProfile("h_diffq22_etaN","",32,-0.8,0.8,tpmin,tpmax,"");
802 h_diffAq22_eta = new TProfile("h_diffAq22_eta","",32,-0.8,0.8,tpmin,tpmax,"");
803 h_diffAq22_etaP = new TProfile("h_diffAq22_etaP","",32,-0.8,0.8,tpmin,tpmax,"");
804 h_diffAq22_etaN = new TProfile("h_diffAq22_etaN","",32,-0.8,0.8,tpmin,tpmax,"");
805 fOutputList->Add(h_diffq22_eta);
806 fOutputList->Add(h_diffq22_etaP);
807 fOutputList->Add(h_diffq22_etaN);
808 fOutputList->Add(h_diffAq22_eta);
809 fOutputList->Add(h_diffAq22_etaP);
810 fOutputList->Add(h_diffAq22_etaN);
811
812 h_diffq32_pt = new TProfile("h_diffq32_pt","",50,0,5,tpmin,tpmax,"");
813 h_diffq32_ptP = new TProfile("h_diffq32_ptP","",50,0,5,tpmin,tpmax,"");
814 h_diffq32_ptN = new TProfile("h_diffq32_ptN","",50,0,5,tpmin,tpmax,"");
815 h_diffAq32_pt = new TProfile("h_diffAq32_pt","",50,0,5,tpmin,tpmax,"");
816 h_diffAq32_ptP = new TProfile("h_diffAq32_ptP","",50,0,5,tpmin,tpmax,"");
817 h_diffAq32_ptN = new TProfile("h_diffAq32_ptN","",50,0,5,tpmin,tpmax,"");
818 fOutputList->Add(h_diffq32_pt);
819 fOutputList->Add(h_diffq32_ptP);
820 fOutputList->Add(h_diffq32_ptN);
821 fOutputList->Add(h_diffAq32_pt);
822 fOutputList->Add(h_diffAq32_ptP);
823 fOutputList->Add(h_diffAq32_ptN);
824 h_diffq32_eta = new TProfile("h_diffq32_eta","",32,-0.8,0.8,tpmin,tpmax,"");
825 h_diffq32_etaP = new TProfile("h_diffq32_etaP","",32,-0.8,0.8,tpmin,tpmax,"");
826 h_diffq32_etaN = new TProfile("h_diffq32_etaN","",32,-0.8,0.8,tpmin,tpmax,"");
827 h_diffAq32_eta = new TProfile("h_diffAq32_eta","",32,-0.8,0.8,tpmin,tpmax,"");
828 h_diffAq32_etaP = new TProfile("h_diffAq32_etaP","",32,-0.8,0.8,tpmin,tpmax,"");
829 h_diffAq32_etaN = new TProfile("h_diffAq32_etaN","",32,-0.8,0.8,tpmin,tpmax,"");
830 fOutputList->Add(h_diffq32_eta);
831 fOutputList->Add(h_diffq32_etaP);
832 fOutputList->Add(h_diffq32_etaN);
833 fOutputList->Add(h_diffAq32_eta);
834 fOutputList->Add(h_diffAq32_etaP);
835 fOutputList->Add(h_diffAq32_etaN);
836
837 h_diffq42_pt = new TProfile("h_diffq42_pt","",50,0,5,tpmin,tpmax,"");
838 h_diffq42_ptP = new TProfile("h_diffq42_ptP","",50,0,5,tpmin,tpmax,"");
839 h_diffq42_ptN = new TProfile("h_diffq42_ptN","",50,0,5,tpmin,tpmax,"");
840 h_diffAq42_pt = new TProfile("h_diffAq42_pt","",50,0,5,tpmin,tpmax,"");
841 h_diffAq42_ptP = new TProfile("h_diffAq42_ptP","",50,0,5,tpmin,tpmax,"");
842 h_diffAq42_ptN = new TProfile("h_diffAq42_ptN","",50,0,5,tpmin,tpmax,"");
843 fOutputList->Add(h_diffq42_pt);
844 fOutputList->Add(h_diffq42_ptP);
845 fOutputList->Add(h_diffq42_ptN);
846 fOutputList->Add(h_diffAq42_pt);
847 fOutputList->Add(h_diffAq42_ptP);
848 fOutputList->Add(h_diffAq42_ptN);
849 h_diffq42_eta = new TProfile("h_diffq42_eta","",32,-0.8,0.8,tpmin,tpmax,"");
850 h_diffq42_etaP = new TProfile("h_diffq42_etaP","",32,-0.8,0.8,tpmin,tpmax,"");
851 h_diffq42_etaN = new TProfile("h_diffq42_etaN","",32,-0.8,0.8,tpmin,tpmax,"");
852 h_diffAq42_eta = new TProfile("h_diffAq42_eta","",32,-0.8,0.8,tpmin,tpmax,"");
853 h_diffAq42_etaP = new TProfile("h_diffAq42_etaP","",32,-0.8,0.8,tpmin,tpmax,"");
854 h_diffAq42_etaN = new TProfile("h_diffAq42_etaN","",32,-0.8,0.8,tpmin,tpmax,"");
855 fOutputList->Add(h_diffq42_eta);
856 fOutputList->Add(h_diffq42_etaP);
857 fOutputList->Add(h_diffq42_etaN);
858 fOutputList->Add(h_diffAq42_eta);
859 fOutputList->Add(h_diffAq42_etaP);
860 fOutputList->Add(h_diffAq42_etaN);
861
862
863
864 h_AT_X_deta = new TProfile("h_AT_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
865 h_AT_X_detaP = new TProfile("h_AT_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
866 h_AT_X_detaN = new TProfile("h_AT_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
867 fOutputList->Add(h_AT_X_deta);
868 fOutputList->Add(h_AT_X_detaP);
869 fOutputList->Add(h_AT_X_detaN);
870 //h_diffq22_X_deta = new TProfile("h_diffq22_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
871 h_diffq22_X_detaP = new TProfile("h_diffq22_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
872 h_diffq22_X_detaN = new TProfile("h_diffq22_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
873 //h_ATdiffq22_X_deta = new TProfile("h_ATdiffq22_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
874 h_ATdiffq22_X_detaP = new TProfile("h_ATdiffq22_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
875 h_ATdiffq22_X_detaN = new TProfile("h_ATdiffq22_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
876 fOutputList->Add(h_diffq22_X_detaP);
877 fOutputList->Add(h_diffq22_X_detaN);
878 fOutputList->Add(h_ATdiffq22_X_detaP);
879 fOutputList->Add(h_ATdiffq22_X_detaN);
880 //h_diffq32_X_deta = new TProfile("h_diffq32_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
881 h_diffq32_X_detaP = new TProfile("h_diffq32_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
882 h_diffq32_X_detaN = new TProfile("h_diffq32_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
883 //h_ATdiffq32_X_deta = new TProfile("h_ATdiffq32_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
884 h_ATdiffq32_X_detaP = new TProfile("h_ATdiffq32_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
885 h_ATdiffq32_X_detaN = new TProfile("h_ATdiffq32_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
886 fOutputList->Add(h_diffq32_X_detaP);
887 fOutputList->Add(h_diffq32_X_detaN);
888 fOutputList->Add(h_ATdiffq32_X_detaP);
889 fOutputList->Add(h_ATdiffq32_X_detaN);
890 //h_diffq42_X_deta = new TProfile("h_diffq42_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
891 h_diffq42_X_detaP = new TProfile("h_diffq42_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
892 h_diffq42_X_detaN = new TProfile("h_diffq42_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
893 //h_ATdiffq42_X_deta = new TProfile("h_ATdiffq42_X_deta","",32,-1.6,1.6,tpmin,tpmax,"");
894 h_ATdiffq42_X_detaP = new TProfile("h_ATdiffq42_X_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
895 h_ATdiffq42_X_detaN = new TProfile("h_ATdiffq42_X_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
896 fOutputList->Add(h_diffq42_X_detaP);
897 fOutputList->Add(h_diffq42_X_detaN);
898 fOutputList->Add(h_ATdiffq42_X_detaP);
899 fOutputList->Add(h_ATdiffq42_X_detaN);
900
901
902 // --- now field selection
903 // --- field 1
904 h_AT_X_deta_F1 = new TProfile("h_AT_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
905 h_AT_X_detaP_F1 = new TProfile("h_AT_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
906 h_AT_X_detaN_F1 = new TProfile("h_AT_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
907 fOutputList->Add(h_AT_X_deta_F1);
908 fOutputList->Add(h_AT_X_detaP_F1);
909 fOutputList->Add(h_AT_X_detaN_F1);
910 //h_diffq22_X_deta_F1 = new TProfile("h_diffq22_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
911 h_diffq22_X_detaP_F1 = new TProfile("h_diffq22_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
912 h_diffq22_X_detaN_F1 = new TProfile("h_diffq22_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
913 //h_ATdiffq22_X_deta_F1 = new TProfile("h_ATdiffq22_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
914 h_ATdiffq22_X_detaP_F1 = new TProfile("h_ATdiffq22_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
915 h_ATdiffq22_X_detaN_F1 = new TProfile("h_ATdiffq22_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
916 fOutputList->Add(h_diffq22_X_detaP_F1);
917 fOutputList->Add(h_diffq22_X_detaN_F1);
918 fOutputList->Add(h_ATdiffq22_X_detaP_F1);
919 fOutputList->Add(h_ATdiffq22_X_detaN_F1);
920 //h_diffq32_X_deta_F1 = new TProfile("h_diffq32_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
921 h_diffq32_X_detaP_F1 = new TProfile("h_diffq32_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
922 h_diffq32_X_detaN_F1 = new TProfile("h_diffq32_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
923 //h_ATdiffq32_X_deta_F1 = new TProfile("h_ATdiffq32_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
924 h_ATdiffq32_X_detaP_F1 = new TProfile("h_ATdiffq32_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
925 h_ATdiffq32_X_detaN_F1 = new TProfile("h_ATdiffq32_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
926 fOutputList->Add(h_diffq32_X_detaP_F1);
927 fOutputList->Add(h_diffq32_X_detaN_F1);
928 fOutputList->Add(h_ATdiffq32_X_detaP_F1);
929 fOutputList->Add(h_ATdiffq32_X_detaN_F1);
930 //h_diffq42_X_deta_F1 = new TProfile("h_diffq42_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
931 h_diffq42_X_detaP_F1 = new TProfile("h_diffq42_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
932 h_diffq42_X_detaN_F1 = new TProfile("h_diffq42_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
933 //h_ATdiffq42_X_deta_F1 = new TProfile("h_ATdiffq42_X_deta_F1","",32,-1.6,1.6,tpmin,tpmax,"");
934 h_ATdiffq42_X_detaP_F1 = new TProfile("h_ATdiffq42_X_detaP_F1","",32,-1.6,1.6,tpmin,tpmax,"");
935 h_ATdiffq42_X_detaN_F1 = new TProfile("h_ATdiffq42_X_detaN_F1","",32,-1.6,1.6,tpmin,tpmax,"");
936 fOutputList->Add(h_diffq42_X_detaP_F1);
937 fOutputList->Add(h_diffq42_X_detaN_F1);
938 fOutputList->Add(h_ATdiffq42_X_detaP_F1);
939 fOutputList->Add(h_ATdiffq42_X_detaN_F1);
940 // --- field 3
941 h_AT_X_deta_F3 = new TProfile("h_AT_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
942 h_AT_X_detaP_F3 = new TProfile("h_AT_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
943 h_AT_X_detaN_F3 = new TProfile("h_AT_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
944 fOutputList->Add(h_AT_X_deta_F3);
945 fOutputList->Add(h_AT_X_detaP_F3);
946 fOutputList->Add(h_AT_X_detaN_F3);
947 //h_diffq22_X_deta_F3 = new TProfile("h_diffq22_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
948 h_diffq22_X_detaP_F3 = new TProfile("h_diffq22_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
949 h_diffq22_X_detaN_F3 = new TProfile("h_diffq22_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
950 //h_ATdiffq22_X_deta_F3 = new TProfile("h_ATdiffq22_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
951 h_ATdiffq22_X_detaP_F3 = new TProfile("h_ATdiffq22_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
952 h_ATdiffq22_X_detaN_F3 = new TProfile("h_ATdiffq22_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
953 fOutputList->Add(h_diffq22_X_detaP_F3);
954 fOutputList->Add(h_diffq22_X_detaN_F3);
955 fOutputList->Add(h_ATdiffq22_X_detaP_F3);
956 fOutputList->Add(h_ATdiffq22_X_detaN_F3);
957 //h_diffq32_X_deta_F3 = new TProfile("h_diffq32_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
958 h_diffq32_X_detaP_F3 = new TProfile("h_diffq32_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
959 h_diffq32_X_detaN_F3 = new TProfile("h_diffq32_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
960 //h_ATdiffq32_X_deta_F3 = new TProfile("h_ATdiffq32_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
961 h_ATdiffq32_X_detaP_F3 = new TProfile("h_ATdiffq32_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
962 h_ATdiffq32_X_detaN_F3 = new TProfile("h_ATdiffq32_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
963 fOutputList->Add(h_diffq32_X_detaP_F3);
964 fOutputList->Add(h_diffq32_X_detaN_F3);
965 fOutputList->Add(h_ATdiffq32_X_detaP_F3);
966 fOutputList->Add(h_ATdiffq32_X_detaN_F3);
967 //h_diffq42_X_deta_F3 = new TProfile("h_diffq42_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
968 h_diffq42_X_detaP_F3 = new TProfile("h_diffq42_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
969 h_diffq42_X_detaN_F3 = new TProfile("h_diffq42_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
970 //h_ATdiffq42_X_deta_F3 = new TProfile("h_ATdiffq42_X_deta_F3","",32,-1.6,1.6,tpmin,tpmax,"");
971 h_ATdiffq42_X_detaP_F3 = new TProfile("h_ATdiffq42_X_detaP_F3","",32,-1.6,1.6,tpmin,tpmax,"");
972 h_ATdiffq42_X_detaN_F3 = new TProfile("h_ATdiffq42_X_detaN_F3","",32,-1.6,1.6,tpmin,tpmax,"");
973 fOutputList->Add(h_diffq42_X_detaP_F3);
974 fOutputList->Add(h_diffq42_X_detaN_F3);
975 fOutputList->Add(h_ATdiffq42_X_detaP_F3);
976 fOutputList->Add(h_ATdiffq42_X_detaN_F3);
977
978
979
980 h_AT_X_dpt = new TProfile("h_AT_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
981 h_AT_X_dptP = new TProfile("h_AT_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
982 h_AT_X_dptN = new TProfile("h_AT_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
983 fOutputList->Add(h_AT_X_dpt);
984 fOutputList->Add(h_AT_X_dptP);
985 fOutputList->Add(h_AT_X_dptN);
986 //h_diffq22_X_dpt = new TProfile("h_diffq22_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
987 h_diffq22_X_dptP = new TProfile("h_diffq22_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
988 h_diffq22_X_dptN = new TProfile("h_diffq22_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
989 //h_ATdiffq22_X_dpt = new TProfile("h_ATdiffq22_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
990 h_ATdiffq22_X_dptP = new TProfile("h_ATdiffq22_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
991 h_ATdiffq22_X_dptN = new TProfile("h_ATdiffq22_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
992 fOutputList->Add(h_diffq22_X_dptP);
993 fOutputList->Add(h_diffq22_X_dptN);
994 fOutputList->Add(h_ATdiffq22_X_dptP);
995 fOutputList->Add(h_ATdiffq22_X_dptN);
996 //h_diffq32_X_dpt = new TProfile("h_diffq32_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
997 h_diffq32_X_dptP = new TProfile("h_diffq32_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
998 h_diffq32_X_dptN = new TProfile("h_diffq32_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
999 //h_ATdiffq32_X_dpt = new TProfile("h_ATdiffq32_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
1000 h_ATdiffq32_X_dptP = new TProfile("h_ATdiffq32_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
1001 h_ATdiffq32_X_dptN = new TProfile("h_ATdiffq32_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
1002 fOutputList->Add(h_diffq32_X_dptP);
1003 fOutputList->Add(h_diffq32_X_dptN);
1004 fOutputList->Add(h_ATdiffq32_X_dptP);
1005 fOutputList->Add(h_ATdiffq32_X_dptN);
1006 //h_diffq42_X_dpt = new TProfile("h_diffq42_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
1007 h_diffq42_X_dptP = new TProfile("h_diffq42_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
1008 h_diffq42_X_dptN = new TProfile("h_diffq42_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
1009 //h_ATdiffq42_X_dpt = new TProfile("h_ATdiffq42_X_dpt","",32,-1.6,1.6,tpmin,tpmax,"");
1010 h_ATdiffq42_X_dptP = new TProfile("h_ATdiffq42_X_dptP","",32,-1.6,1.6,tpmin,tpmax,"");
1011 h_ATdiffq42_X_dptN = new TProfile("h_ATdiffq42_X_dptN","",32,-1.6,1.6,tpmin,tpmax,"");
1012 fOutputList->Add(h_diffq42_X_dptP);
1013 fOutputList->Add(h_diffq42_X_dptN);
1014 fOutputList->Add(h_ATdiffq42_X_dptP);
1015 fOutputList->Add(h_ATdiffq42_X_dptN);
1016
1017
1018
1019 // --- now Xsub...
1020
1021 h_AT_S_deta = new TProfile("h_AT_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1022 h_AT_S_detaP = new TProfile("h_AT_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1023 h_AT_S_detaN = new TProfile("h_AT_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1024 fOutputList->Add(h_AT_S_deta);
1025 fOutputList->Add(h_AT_S_detaP);
1026 fOutputList->Add(h_AT_S_detaN);
1027 //h_diffq22_S_deta = new TProfile("h_diffq22_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1028 h_diffq22_S_detaP = new TProfile("h_diffq22_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1029 h_diffq22_S_detaN = new TProfile("h_diffq22_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1030 //h_ATdiffq22_S_deta = new TProfile("h_ATdiffq22_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1031 h_ATdiffq22_S_detaP = new TProfile("h_ATdiffq22_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1032 h_ATdiffq22_S_detaN = new TProfile("h_ATdiffq22_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1033 fOutputList->Add(h_diffq22_S_detaP);
1034 fOutputList->Add(h_diffq22_S_detaN);
1035 fOutputList->Add(h_ATdiffq22_S_detaP);
1036 fOutputList->Add(h_ATdiffq22_S_detaN);
1037 //h_diffq32_S_deta = new TProfile("h_diffq32_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1038 h_diffq32_S_detaP = new TProfile("h_diffq32_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1039 h_diffq32_S_detaN = new TProfile("h_diffq32_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1040 //h_ATdiffq32_S_deta = new TProfile("h_ATdiffq32_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1041 h_ATdiffq32_S_detaP = new TProfile("h_ATdiffq32_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1042 h_ATdiffq32_S_detaN = new TProfile("h_ATdiffq32_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1043 fOutputList->Add(h_diffq32_S_detaP);
1044 fOutputList->Add(h_diffq32_S_detaN);
1045 fOutputList->Add(h_ATdiffq32_S_detaP);
1046 fOutputList->Add(h_ATdiffq32_S_detaN);
1047 //h_diffq42_S_deta = new TProfile("h_diffq42_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1048 h_diffq42_S_detaP = new TProfile("h_diffq42_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1049 h_diffq42_S_detaN = new TProfile("h_diffq42_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1050 //h_ATdiffq42_S_deta = new TProfile("h_ATdiffq42_S_deta","",32,-1.6,1.6,tpmin,tpmax,"");
1051 h_ATdiffq42_S_detaP = new TProfile("h_ATdiffq42_S_detaP","",32,-1.6,1.6,tpmin,tpmax,"");
1052 h_ATdiffq42_S_detaN = new TProfile("h_ATdiffq42_S_detaN","",32,-1.6,1.6,tpmin,tpmax,"");
1053 fOutputList->Add(h_diffq42_S_detaP);
1054 fOutputList->Add(h_diffq42_S_detaN);
1055 fOutputList->Add(h_ATdiffq42_S_detaP);
1056 fOutputList->Add(h_ATdiffq42_S_detaN);
1057
1058 // ---
1059
1060 // --- random fun
1061
1062 // ---
1063
1064 h_a1q22_cent = new TProfile("h_a1q22_cent","q22 vs centrality",100,0,100,tpmin,tpmax,"");
1065 h_a1q22_centP = new TProfile("h_a1q22_centP","q22(pos) vs centrality",100,0,100,tpmin,tpmax,"");
1066 h_a1q22_centN = new TProfile("h_a1q22_centN","q22(neg) vs centrality",100,0,100,tpmin,tpmax,"");
1067 fOutputList->Add(h_a1q22_cent);
1068 fOutputList->Add(h_a1q22_centP);
1069 fOutputList->Add(h_a1q22_centN);
1070
1071 h_a2q22_cent = new TProfile("h_a2q22_cent","q22 vs centrality",100,0,100,tpmin,tpmax,"");
1072 h_a2q22_centP = new TProfile("h_a2q22_centP","q22(pos) vs centrality",100,0,100,tpmin,tpmax,"");
1073 h_a2q22_centN = new TProfile("h_a2q22_centN","q22(neg) vs centrality",100,0,100,tpmin,tpmax,"");
1074 fOutputList->Add(h_a2q22_cent);
1075 fOutputList->Add(h_a2q22_centP);
1076 fOutputList->Add(h_a2q22_centN);
1077
1078 h_rAq22_X1_cent = new TProfile("h_rAq22_X1_cent","q22 vs centrality",100,0,100,tpmin,tpmax,"");
1079 h_rAq22_X1_centP = new TProfile("h_rAq22_X1_centP","q22(pos) vs centrality",100,0,100,tpmin,tpmax,"");
1080 h_rAq22_X1_centN = new TProfile("h_rAq22_X1_centN","q22(neg) vs centrality",100,0,100,tpmin,tpmax,"");
1081 h_rAq22_X2_cent = new TProfile("h_rAq22_X2_cent","q22 vs centrality",100,0,100,tpmin,tpmax,"");
1082 h_rAq22_X2_centP = new TProfile("h_rAq22_X2_centP","q22(pos) vs centrality",100,0,100,tpmin,tpmax,"");
1083 h_rAq22_X2_centN = new TProfile("h_rAq22_X2_centN","q22(neg) vs centrality",100,0,100,tpmin,tpmax,"");
1084 h_rAq22_X3_cent = new TProfile("h_rAq22_X3_cent","q22 vs centrality",100,0,100,tpmin,tpmax,"");
1085 h_rAq22_X3_centP = new TProfile("h_rAq22_X3_centP","q22(pos) vs centrality",100,0,100,tpmin,tpmax,"");
1086 h_rAq22_X3_centN = new TProfile("h_rAq22_X3_centN","q22(neg) vs centrality",100,0,100,tpmin,tpmax,"");
1087 h_rAq22_X4_cent = new TProfile("h_rAq22_X4_cent","q22 vs centrality",100,0,100,tpmin,tpmax,"");
1088 h_rAq22_X4_centP = new TProfile("h_rAq22_X4_centP","q22(pos) vs centrality",100,0,100,tpmin,tpmax,"");
1089 h_rAq22_X4_centN = new TProfile("h_rAq22_X4_centN","q22(neg) vs centrality",100,0,100,tpmin,tpmax,"");
1090 fOutputList->Add(h_rAq22_X1_cent);
1091 fOutputList->Add(h_rAq22_X1_centP);
1092 fOutputList->Add(h_rAq22_X1_centN);
1093 fOutputList->Add(h_rAq22_X2_cent);
1094 fOutputList->Add(h_rAq22_X2_centP);
1095 fOutputList->Add(h_rAq22_X2_centN);
1096 fOutputList->Add(h_rAq22_X3_cent);
1097 fOutputList->Add(h_rAq22_X3_centP);
1098 fOutputList->Add(h_rAq22_X3_centN);
1099 fOutputList->Add(h_rAq22_X4_cent);
1100 fOutputList->Add(h_rAq22_X4_centP);
1101 fOutputList->Add(h_rAq22_X4_centN);
1102
1103
1104
1105 // ---------------------------- //
1106 // --- done with histograms --- //
1107 // ---------------------------- //
1108
1109}
1110
1111
1112
1113//----------------------------------------------
1114void AliAnalysisTaskCMEv2A::UserExec(Option_t *)
1115{
1116
1117 // Main analyis loop, called for each event
1118
1119 if(debug>0) cout<<"Processing event with debug = "<<debug<<endl;
1120
1121 // Analysis manager
1122 AliAnalysisManager *manager = AliAnalysisManager::GetAnalysisManager();
1123 if(!manager)
1124 {
1125 cout<<"FATAL: could not get Analysis Manager."<<endl;
1126 return;
1127 }
1128
1129 // Input handler
1130 AliAODInputHandler *handler = (AliAODInputHandler *)manager->GetInputEventHandler();
1131 if(!handler)
1132 {
1133 cout<<"FATAL: could not get Input Handler."<<endl;
1134 return;
1135 }
1136
1137
1138 // AOD Event object from tree
1139 AliAODEvent *fAOD = (AliAODEvent *)InputEvent();
1140 if(!fAOD)
1141 {
1142 if(debug>-1) cout<<"ERROR: AOD event object not available. Discarding event..."<<endl;
1143 return;
1144 }
1145
1146
1147 int runnumber = fAOD->GetRunNumber();
1148 float mag = fAOD->GetMagneticField();
1149
1150 if(debug>0)
1151 {
1152 cout<<"runnumber is "<<runnumber<<endl;
1153 cout<<"magnetic field is "<<mag<<endl;
1154 }
1155
1156 // MC Event object from tree
1157 AliMCEvent *fMC = MCEvent();
1158 if(!fMC&&doMC)
1159 {
1160 if(debug>-1) cout<<"ERROR: MC event object not available. Discarding event..."<<endl;
1161 return;
1162 }
1163
1164 // get pid
1165 AliPIDResponse *fPID = handler->GetPIDResponse();
2d95bc0f 1166 if(!fPID && !doMC) // use PID only if no MC
545d70ee 1167 {
1168 if(debug>-1) cout<<"ERROR: PIDResponse object not available. Discarding event..."<<endl;
1169 return;
1170 }
1171
1172 // // Eventplane object (from AOD)
1173 // AliEventplane *fEventplane = fAOD->GetEventplane();
1174 // if(!fEventplane)
1175 // {
1176 // if(debug>-1) cout<<"ERROR: Eventplane object not available. Discarding event..."<<endl;
1177 // return;
1178 // }
1179
1180 // float psi_V0_h2 = fEventplane->GetEventplane("V0",fAOD,2);
1181 // float psi_V0A_h2 = fEventplane->GetEventplane("V0A",fAOD,2);
1182 // float psi_V0C_h2 = fEventplane->GetEventplane("V0C",fAOD,2);
1183
1184 // fHistPlaneV0h2->Fill(psi_V0_h2);
1185 // fHistPlaneV0Ah2->Fill(psi_V0A_h2);
1186 // fHistPlaneV0Ch2->Fill(psi_V0C_h2);
1187 // fHistPlaneV0ACDCh2->Fill(psi_V0A_h2-psi_V0C_h2);
1188
1189
1190
1191 // Centrality object (from AOD)
1192 AliCentrality *fCentrality = fAOD->GetCentrality();
1193 if(!fCentrality)
1194 {
1195 if(debug>-1) cout<<"ERROR: Centrality object not available. Discarding event..."<<endl;
1196 return;
1197 }
1198
1199 float centTRK = fCentrality->GetCentralityPercentile("TRK");
1200 float centV0M = fCentrality->GetCentralityPercentile("V0M");
1201 float centSPD = fCentrality->GetCentralityPercentile("CL1");//outer SPD?
1202 float cent = centTRK;
1203 if(centhandle==2) cent = centV0M;
1204 if(centhandle==3) cent = centSPD;
1205
1206 int icent = int(cent)/10;
1207
1208 int centstatus = 0;
1209 if(centTRK<0.0||centV0M<0.0) centstatus = -1;
1210 if(centTRK<0.0&&centV0M<0.0) centstatus = -2;
1211 if(centTRK>0.0||centV0M>0.0) centstatus = 1;
1212 if(centTRK>0.0&&centV0M>0.0) centstatus = 2;
1213 if(centTRK==0.0&&centV0M==0.0) centstatus = 0;
1214
1215 if(debug>0)
1216 {
1217 cout<<"centTRK "<<centTRK<<endl;
1218 cout<<"centV0M "<<centV0M<<endl;
1219 cout<<"centSPD "<<centSPD<<endl;
1220 cout<<"centrality selection is "<<centhandle<<endl;
1221 cout<<"cent is "<<cent<<endl;
1222 }
1223 fHistCentTRK->Fill(centTRK);
1224 fHistCentV0M->Fill(centV0M);
1225 fHistCentDIFF->Fill(centTRK-centV0M);
1226
1227 //ULong64_t mask = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->IsEventSelected();
1228 ULong64_t mask = handler->IsEventSelected();
1229 ULong64_t amb = AliVEvent::kMB;
1230 ULong64_t acn = AliVEvent::kCentral;
1231 ULong64_t asc = AliVEvent::kSemiCentral;
1232
1233 if(debug>0) cout<<"trigger selection is "<<trigger<<endl;
1234 if(debug>0) cout<<"trigger mask is "<<mask<<endl;
1235
1236 if(mask&amb)
1237 {
1238 fHistCentTRKAVEkMB->Fill(centTRK);
1239 fHistCentV0MAVEkMB->Fill(centV0M);
1240 }
1241 if(mask&acn)
1242 {
1243 fHistCentTRKAVEkCentral->Fill(centTRK);
1244 fHistCentV0MAVEkCentral->Fill(centV0M);
1245 }
1246 if(mask&asc)
1247 {
1248 fHistCentTRKAVEkSemiCentral->Fill(centTRK);
1249 fHistCentV0MAVEkSemiCentral->Fill(centV0M);
1250 }
1251 if(mask&(amb|acn|asc))
1252 {
1253 fHistCentTRKAVEkA3->Fill(centTRK);
1254 fHistCentV0MAVEkA3->Fill(centV0M);
1255 }
1256
1257 if(mask&trigger)
1258 {
1259 fHistCentTRKAVEkSel->Fill(centTRK);
1260 fHistCentV0MAVEkSel->Fill(centV0M);
1261 }
1262 else
1263 {
1264 if(debug>0) cout<<"wrong trigger, rejecting event"<<endl;
1265 return;
1266 }
1267
1268 if(fabs(centTRK-centV0M)>centcut)
1269 {
1270 if(debug>0) cout<<"centrality difference outside cut, rejecting event"<<endl;
1271 return;
1272 }
1273
1274
1275
1276 // AOD vertex objects
1277 AliAODVertex *fVtx = fAOD->GetPrimaryVertex();
1278 AliAODVertex *fVtxSPD = fAOD->GetPrimaryVertexSPD();
1279 if(!fVtx)
1280 {
1281 if(debug>-1) cout<<"ERROR: Vertex object not available. Discarding event..."<<endl;
1282 return;
1283 }
1284 float zvtxV0 = fVtx->GetZ();
1285 float zvtxSPD = fVtxSPD->GetZ();
1286 if(debug>0) cout<<"zvtxV0 is "<<zvtxV0<<endl;
1287 if(debug>0) cout<<"zvtxSPD is "<<zvtxSPD<<endl;
1288 if(centstatus==2)
1289 {
1290 fHistZVtx->Fill(zvtxV0);
1291 fHistZVtxD->Fill(zvtxV0-zvtxSPD);
1292 }
1293 if(fabs(zvtxV0)>zvtxcut)
1294 {
1295 if(debug>0) cout<<"vertex outside cut, rejecting event"<<endl;
1296 return;
1297 }
1298 float eventX = fVtx->GetX();
1299 float eventY = fVtx->GetY();
1300 float eventZ = fVtx->GetZ();
1301
1302
1303
1304 // // get V0 object
1305 // AliAODVZERO *fV0 = fAOD->GetVZEROData();
1306 // for(int i=0; i<64; i++)
1307 // {
1308 // float phiV0 = pi/4.0*(0.5+i%8);
1309 // float multV0 = fV0->GetMultiplicity(i);
1310 // }
1311
1312
1313 int d_ntrk = fAOD->GetNumberOfTracks();
1314 int d_ntrkMC = 0;
1315 if(fMC) d_ntrkMC = fMC->GetNumberOfTracks();
1316
1317 if(centstatus==-2&&d_ntrk>0) centstatus = -3;
1318 if(debug>0) cout<<"there are "<<d_ntrk<<" tracks in this event"<<endl;
1319 if(debug>0&&fMC) cout<<"there are "<<d_ntrkMC<<" Monte Carlo tracks in this event"<<endl;
1320 if(debug>0) cout<<"centrality diagnostic is "<<centstatus<<endl;
1321
1322 fHistCentDIAG->Fill(centstatus);
1323 if(centstatus==-3)
1324 {
1325 fHistCtrkDIAG->Fill(d_ntrk);
1326 fHistVtxRDIAG->Fill(zvtxV0);
1327 fHistVtxSDIAG->Fill(zvtxSPD);
1328 }
1329 if(centstatus==-2)
1330 {
1331 fHistVtxRDIBG->Fill(zvtxV0);
1332 fHistVtxSDIBG->Fill(zvtxSPD);
1333 }
1334
1335 // --- AliAnalysisUtils for pileup cuts
1336 if(dopupcut)
1337 {
1338 AliAnalysisUtils *fUtils = new AliAnalysisUtils();
1339 if(!fUtils)
1340 {
1341 if(debug>-1) cout<<"ERROR: cannot find AliAnalysisUtils..."<<endl;
1342 return;
1343 }
1344 fUtils->SetUseOutOfBunchPileUp(true);
1345 bool pileup = fUtils->IsPileUpEvent(fAOD);
1346 //bool pileup = fAOD->IsPileUpFromSPD(3,0.8,3.0,2.0,5.0); // default parameters in AliAODEvent.h
1347 if(pileup)
1348 {
1349 if(debug>0) cout<<"Rejecting event for pileup (AliAnalysisUtils)"<<endl;
1350 return;
1351 }
1352 if(fUtils) delete fUtils;
1353 //
1354 // ---
1355 //
1356 /*
1357 const AliAODVertex* vtxSPD = fAOD->GetPrimaryVertexSPD();
1358 if(!vtxSPD||vtxSPD->GetNContributors()<=0)
1359 {
1360 if(debug>0) cout<<"rejecting pileup event (no vtxSPD or zero contributors) "<<vtxSPD<<endl;
1361 return;
1362 }
1363
1364 const AliAODVertex* vtxTPC = 0;
1365 int nVertices = fAOD->GetNumberOfVertices();
1366 if(debug>0) cout<<"number of vertices is "<<nVertices<<endl;
1367 for(int iVertices = 0; iVertices < nVertices; iVertices++)
1368 {
1369 const AliAODVertex* vertex = fAOD->GetVertex(iVertices);
1370 if(debug>0) cout<<"vertex type is "<<vertex->GetType()<<endl;
1371 if(vertex->GetType()!=AliAODVertex::kMainTPC) continue;
1372 vtxTPC = vertex;
1373 }
1374 if(!vtxTPC||vtxTPC->GetNContributors()<=0)
1375 {
1376 if(debug>0) cout<<"rejecting pileup event (no vtxTPC or zero contributors)"<<vtxTPC<<endl;
1377 return;
1378 }
1379
1380 float diffZ = vtxSPD->GetZ() - vtxTPC->GetZ();
1381 if(fabs(diffZ)>2.0)
1382 {
1383 if(debug>0) cout<<"rejecting pileup event with vtxTPC "<<vtxTPC->GetZ()<<" vtxSPD "<<vtxSPD->GetZ()<<endl;
1384 return;
1385 }
1386 */
1387 //
1388 // ---
1389 //
1390 }
1391
1392
1393 int ntrk = 0;
1394 int ntrkpos = 0;
1395 int ntrkneg = 0;
1396 int ntrkL = 0;
1397 int ntrkposL = 0;
1398 int ntrknegL = 0;
1399 int ntrkR = 0;
1400 int ntrkposR = 0;
1401 int ntrknegR = 0;
1402 int ntrkA1 = 0;
1403 int ntrkposA1 = 0;
1404 int ntrknegA1 = 0;
1405 int ntrkA2 = 0;
1406 int ntrkposA2 = 0;
1407 int ntrknegA2 = 0;
1408 int ntrkMC = 0;
1409 int ntrkposMC = 0;
1410 int ntrknegMC = 0;
1411
1412 int cutntrk = 0;
1413 int cutntrkpos = 0;
1414 int cutntrkneg = 0;
1415 int cutntrkL = 0;
1416 int cutntrkposL = 0;
1417 int cutntrknegL = 0;
1418 int cutntrkR = 0;
1419 int cutntrkposR = 0;
1420 int cutntrknegR = 0;
1421 int cutntrkMC = 0;
1422 int cutntrkposMC = 0;
1423 int cutntrknegMC = 0;
1424
1425 const int hmax = 9; // number of harmonics including 0 for multiplicity counting
1426 float tpcXplo[hmax], tpcYplo[hmax], tpcXpro[hmax], tpcYpro[hmax];
1427 float tpcXnlo[hmax], tpcYnlo[hmax], tpcXnro[hmax], tpcYnro[hmax];
1428 float tpcXpli[hmax], tpcYpli[hmax], tpcXpri[hmax], tpcYpri[hmax];
1429 float tpcXnli[hmax], tpcYnli[hmax], tpcXnri[hmax], tpcYnri[hmax];
1430 float tpcXpl[hmax], tpcYpl[hmax], tpcXpr[hmax], tpcYpr[hmax];
1431 float tpcXnl[hmax], tpcYnl[hmax], tpcXnr[hmax], tpcYnr[hmax];
1432 for(int i=0; i<hmax;i++)
1433 {
1434 tpcXplo[i] = 0.0;
1435 tpcYplo[i] = 0.0;
1436 tpcXpro[i] = 0.0;
1437 tpcYpro[i] = 0.0;
1438 tpcXnlo[i] = 0.0;
1439 tpcYnlo[i] = 0.0;
1440 tpcXnro[i] = 0.0;
1441 tpcYnro[i] = 0.0;
1442 //
1443 tpcXpli[i] = 0.0;
1444 tpcYpli[i] = 0.0;
1445 tpcXpri[i] = 0.0;
1446 tpcYpri[i] = 0.0;
1447 tpcXnli[i] = 0.0;
1448 tpcYnli[i] = 0.0;
1449 tpcXnri[i] = 0.0;
1450 tpcYnri[i] = 0.0;
1451 //
1452 tpcXpl[i] = 0.0;
1453 tpcYpl[i] = 0.0;
1454 tpcXpr[i] = 0.0;
1455 tpcYpr[i] = 0.0;
1456 tpcXnl[i] = 0.0;
1457 tpcYnl[i] = 0.0;
1458 tpcXnr[i] = 0.0;
1459 tpcYnr[i] = 0.0;
1460 }
1461
1462 float MCtpcXplo[hmax], MCtpcYplo[hmax], MCtpcXpro[hmax], MCtpcYpro[hmax];
1463 float MCtpcXnlo[hmax], MCtpcYnlo[hmax], MCtpcXnro[hmax], MCtpcYnro[hmax];
1464 float MCtpcXpli[hmax], MCtpcYpli[hmax], MCtpcXpri[hmax], MCtpcYpri[hmax];
1465 float MCtpcXnli[hmax], MCtpcYnli[hmax], MCtpcXnri[hmax], MCtpcYnri[hmax];
1466 float MCtpcXpl[hmax], MCtpcYpl[hmax], MCtpcXpr[hmax], MCtpcYpr[hmax];
1467 float MCtpcXnl[hmax], MCtpcYnl[hmax], MCtpcXnr[hmax], MCtpcYnr[hmax];
1468 for(int i=0; i<hmax;i++)
1469 {
1470 MCtpcXplo[i] = 0.0;
1471 MCtpcYplo[i] = 0.0;
1472 MCtpcXpro[i] = 0.0;
1473 MCtpcYpro[i] = 0.0;
1474 MCtpcXnlo[i] = 0.0;
1475 MCtpcYnlo[i] = 0.0;
1476 MCtpcXnro[i] = 0.0;
1477 MCtpcYnro[i] = 0.0;
1478 //
1479 MCtpcXpli[i] = 0.0;
1480 MCtpcYpli[i] = 0.0;
1481 MCtpcXpri[i] = 0.0;
1482 MCtpcYpri[i] = 0.0;
1483 MCtpcXnli[i] = 0.0;
1484 MCtpcYnli[i] = 0.0;
1485 MCtpcXnri[i] = 0.0;
1486 MCtpcYnri[i] = 0.0;
1487 //
1488 MCtpcXpl[i] = 0.0;
1489 MCtpcYpl[i] = 0.0;
1490 MCtpcXpr[i] = 0.0;
1491 MCtpcYpr[i] = 0.0;
1492 MCtpcXnl[i] = 0.0;
1493 MCtpcYnl[i] = 0.0;
1494 MCtpcXnr[i] = 0.0;
1495 MCtpcYnr[i] = 0.0;
1496 }
1497
1498 // ---------------------------------- //
1499 // --- Now looping over MC tracks --- //
1500 // ---------------------------------- //
1501 if(fMC)
1502 {
1503 // variables for 3rd particle
1504 float MCpt3[d_ntrkMC];
1505 float MCeta3[d_ntrkMC];
1506 //float MCphi3[d_ntrkMC];//new...
1507 int MCcharge3[d_ntrkMC];
1508 // initial variables to out of range values
1509 // for cases when the track loop skips certain tracks
1510 for(int i=0; i<d_ntrk; i++)
1511 {
1512 MCpt3[i] = -99;
1513 MCeta3[i] = -99;
1514 //MCphi3[i] = -99;
1515 MCcharge3[i] = 0;
1516 }
1517 // --- first MC track loop
1518 for(int itrkMC = 0; itrkMC<d_ntrkMC; itrkMC++)
1519 {
1520 //AliVParticle *trackMC = fMC->GetTrack(itrkMC);
1521 AliAODMCParticle *trackMC = (AliAODMCParticle *)fMC->GetTrack(itrkMC);
1522 if(!trackMC)
1523 {
1524 if(debug>0) cout<<"ERROR: Could not receive track "<<itrkMC<<" (mc loop)"<<endl;
1525 continue;
1526 }
1527
1528 float pt = trackMC->Pt();
1529 fHistPtMC->Fill(pt);
1530
1531 float phi = trackMC->Phi();
1532 float eta = trackMC->Eta();
1533 int charge = trackMC->Charge();
1534 if(charge==0) continue; // new
1535 charge /= 3; // new
1536 // int ncls = trackMC->GetMCTPCNcls();
1537 // float dedx = trackMC->GetMCTPCsignal();
1538 // float dcaxy = trackMC->DCA();
1539 // float dcaz = trackMC->ZAtDCA();
1540 bool pos = charge>0.0;
1541 bool neg = charge<0.0;
1542
1543 fHistPhiMC->Fill(phi);
1544 fHistEtaMC->Fill(eta);
1545 fHistChargeMC->Fill(charge);
1546 // fHistMCTPCnclsMC->Fill(ncls);
1547 // fHistDedxMC->Fill(dedx);
1548 // fHistDCAxyMC->Fill(dcaxy);
1549 // fHistDCAzMC->Fill(dcaz);
1550
1551 // apply kinematic cuts
1552 if(fabs(eta)>outeta) continue;
1553 if(pt<ptmin||pt>ptmax) continue;
1554
1555 //if(charge!=-3&&charge!=+3) continue;// x3 by convention
1556 if(!trackMC->IsPrimary()) continue;
1557 if(!trackMC->IsPhysicalPrimary()) continue;
1558 //if(abs(p0->GetPdgCode())==11) continue; //electrons
1559 //if(abs(p0->GetPdgCode())==13) continue; //electrons
1560
1561 ntrkMC++;
1562 if(pos) ntrkposMC++;
1563 if(neg) ntrknegMC++;
1564
1565 if(pt>ptmin&&pt<ptmax&&fabs(eta)<outeta&&fabs(eta)>excleta)
1566 {
1567 cutntrkMC++;
1568 if(pos) cutntrkposMC++;
1569 if(neg) cutntrknegMC++;
1570 }
1571
1572 if(pt>ptmin&&pt<ptmax)
1573 {
1574 // --- outer
1575 if(eta>-outeta&&eta<-ineta)
1576 {
1577 for(int i=1; i<hmax;i++)
1578 {
1579 if(charge>0)
1580 {
1581 MCtpcXplo[i] += cos(i*phi);
1582 MCtpcYplo[i] += sin(i*phi);
1583 MCtpcXpl[i] += cos(i*phi);
1584 MCtpcYpl[i] += sin(i*phi);
1585 }
1586 if(charge<0)
1587 {
1588 MCtpcXnlo[i] += cos(i*phi);
1589 MCtpcYnlo[i] += sin(i*phi);
1590 MCtpcXnl[i] += cos(i*phi);
1591 MCtpcYnl[i] += sin(i*phi);
1592 }
1593 }
1594 if(charge>0)
1595 {
1596 MCtpcXplo[0] += 1.0;
1597 MCtpcYplo[0] += pt;
1598 MCtpcXpl[0] += 1.0;
1599 MCtpcYpl[0] += pt;
1600 }
1601 if(charge<0)
1602 {
1603 MCtpcXnlo[0] += 1.0;
1604 MCtpcYnlo[0] += pt;
1605 MCtpcXnl[0] += 1.0;
1606 MCtpcYnl[0] += pt;
1607 }
1608 } // end left half of MCTPC
1609 //----------------------------
1610 if(eta>ineta&&eta<outeta)
1611 {
1612 for(int i=1; i<hmax; i++)
1613 {
1614 if(charge>0)
1615 {
1616 MCtpcXpro[i] += cos(i*phi);
1617 MCtpcYpro[i] += sin(i*phi);
1618 MCtpcXpr[i] += cos(i*phi);
1619 MCtpcYpr[i] += sin(i*phi);
1620 }
1621 if(charge<0)
1622 {
1623 MCtpcXnro[i] += cos(i*phi);
1624 MCtpcYnro[i] += sin(i*phi);
1625 MCtpcXnr[i] += cos(i*phi);
1626 MCtpcYnr[i] += sin(i*phi);
1627 }
1628 }
1629 if(charge>0)
1630 {
1631 MCtpcXpro[0] += 1.0;
1632 MCtpcYpro[0] += pt;
1633 MCtpcXpr[0] += 1.0;
1634 MCtpcYpr[0] += pt;
1635 }
1636 if(charge<0)
1637 {
1638 MCtpcXnro[0] += 1.0;
1639 MCtpcYnro[0] += pt;
1640 MCtpcXnr[0] += 1.0;
1641 MCtpcYnr[0] += pt;
1642 }
1643 } // end right half of MCTPC
1644 // --- inner
1645 if(eta>-ineta&&eta<-excleta)
1646 {
1647 for(int i=1; i<hmax;i++)
1648 {
1649 if(charge>0)
1650 {
1651 MCtpcXpli[i] += cos(i*phi);
1652 MCtpcYpli[i] += sin(i*phi);
1653 MCtpcXpl[i] += cos(i*phi);
1654 MCtpcYpl[i] += sin(i*phi);
1655 }
1656 if(charge<0)
1657 {
1658 MCtpcXnli[i] += cos(i*phi);
1659 MCtpcYnli[i] += sin(i*phi);
1660 MCtpcXnl[i] += cos(i*phi);
1661 MCtpcYnl[i] += sin(i*phi);
1662 }
1663 }
1664 if(charge>0)
1665 {
1666 MCtpcXpli[0] += 1.0;
1667 MCtpcYpli[0] += pt;
1668 MCtpcXpl[0] += 1.0;
1669 MCtpcYpl[0] += pt;
1670 }
1671 if(charge<0)
1672 {
1673 MCtpcXnli[0] += 1.0;
1674 MCtpcYnli[0] += pt;
1675 MCtpcXnl[0] += 1.0;
1676 MCtpcYnl[0] += pt;
1677 }
1678 } // end left half of MCTPC
1679 //----------------------------
1680 if(eta>excleta&&eta<ineta)
1681 {
1682 for(int i=1; i<hmax; i++)
1683 {
1684 if(charge>0)
1685 {
1686 MCtpcXpri[i] += cos(i*phi);
1687 MCtpcYpri[i] += sin(i*phi);
1688 MCtpcXpr[i] += cos(i*phi);
1689 MCtpcYpr[i] += sin(i*phi);
1690 }
1691 if(charge<0)
1692 {
1693 MCtpcXnri[i] += cos(i*phi);
1694 MCtpcYnri[i] += sin(i*phi);
1695 MCtpcXnr[i] += cos(i*phi);
1696 MCtpcYnr[i] += sin(i*phi);
1697 }
1698 }
1699 if(charge>0)
1700 {
1701 MCtpcXpri[0] += 1.0;
1702 MCtpcYpri[0] += pt;
1703 MCtpcXpr[0] += 1.0;
1704 MCtpcYpr[0] += pt;
1705 }
1706 if(charge<0)
1707 {
1708 MCtpcXnri[0] += 1.0;
1709 MCtpcYnri[0] += pt;
1710 MCtpcXnr[0] += 1.0;
1711 MCtpcYnr[0] += pt;
1712 }
1713 } // end right half of MCTPC
1714 } // end pt selection
1715
1716 h_AT_etaMC->Fill(cent,charge);
1717 h2_AT_etaMC->Fill(cent,charge);
1718
1719 MCpt3[itrkMC] = pt;
1720 MCeta3[itrkMC] = eta;
1721 //MCphi3[itrkMC] = phi;
1722 MCcharge3[itrkMC] = charge;
1723
1724 } // end of first MC track loop
1725
1726 // if(debug>0)
1727 // {
1728 // cout<<"number of tracks with wrong filter bit "<<badbit<<endl;
1729 // cout<<"difference = ntrk - nbad = "<<d_ntrk-badbit<<endl;
1730 // }
1731 float MCtpcX[9], MCtpcY[9], MCtpcQQ[9];//, qq[9];
1732 float MCtpcXp[9], MCtpcYp[9], MCtpcQQp[9];//, qqp[9]; // pos
1733 float MCtpcXn[9], MCtpcYn[9], MCtpcQQn[9];//, qqn[9]; // neg
1734
1735 float MCqasymm = float(ntrkposMC-ntrknegMC)/float(ntrkMC);
1736 if(doacuts)
1737 {
1738 MCqasymm = float(cutntrkposMC-cutntrknegMC)/float(cutntrkMC);
1739 }
1740 h_A_centMC->Fill(cent,MCqasymm);
1741 h2_A_centMC->Fill(cent,MCqasymm);
1742
1743
1744 for(int i=0; i<6; i++)
1745 {
1746 MCtpcX[i]=MCtpcXpl[i]+MCtpcXnl[i]+MCtpcXpr[i]+MCtpcXnr[i];
1747 MCtpcY[i]=MCtpcYpl[i]+MCtpcYnl[i]+MCtpcYpr[i]+MCtpcYnr[i];
1748 MCtpcQQ[i]=MCtpcX[i]*MCtpcX[i]+MCtpcY[i]*MCtpcY[i];
1749 //qq[i]=sqrt(MCtpcQQ[i]/MCtpcX[0]);
1750 // pos
1751 MCtpcXp[i]=MCtpcXpl[i]+MCtpcXpr[i];
1752 MCtpcYp[i]=MCtpcYpl[i]+MCtpcYpr[i];
1753 MCtpcQQp[i]=MCtpcXp[i]*MCtpcXp[i]+MCtpcYp[i]*MCtpcYp[i];
1754 //qqp[i]=sqrt(MCtpcQQp[i]/MCtpcXp[0]);
1755 // neg
1756 MCtpcXn[i]=MCtpcXnl[i]+MCtpcXnr[i];
1757 MCtpcYn[i]=MCtpcYnl[i]+MCtpcYnr[i];
1758 MCtpcQQn[i]=MCtpcXn[i]*MCtpcXn[i]+MCtpcYn[i]*MCtpcYn[i];
1759 //qqn[i]=sqrt(MCtpcQQn[i]/MCtpcXn[0]);
1760 }
1761
1762 float M = MCtpcX[0];
1763 float W_2 = M*(M-1);
1764 float Mp = MCtpcXp[0];
1765 float Wp_2 = Mp*(Mp-1);
1766 float Mn = MCtpcXn[0];
1767 float Wn_2 = Mn*(Mn-1);
1768
1769 float MCtpcXl2 = MCtpcXnl[2]+MCtpcXpl[2];
1770 float MCtpcYl2 = MCtpcYnl[2]+MCtpcYpl[2];
1771 float MCtpcXr2 = MCtpcXnr[2]+MCtpcXpr[2];
1772 float MCtpcYr2 = MCtpcYnr[2]+MCtpcYpr[2];
1773 float MCtpcXl0 = MCtpcXnl[0]+MCtpcXpl[0];
1774 float MCtpcXr0 = MCtpcXnr[0]+MCtpcXpr[0];
1775
1776 float MCtpcXl2o = MCtpcXnlo[2]+MCtpcXplo[2];
1777 float MCtpcYl2o = MCtpcYnlo[2]+MCtpcYplo[2];
1778 float MCtpcXr2o = MCtpcXnro[2]+MCtpcXpro[2];
1779 float MCtpcYr2o = MCtpcYnro[2]+MCtpcYpro[2];
1780 float MCtpcXl0o = MCtpcXnlo[0]+MCtpcXplo[0];
1781 float MCtpcXr0o = MCtpcXnro[0]+MCtpcXpro[0];
1782
1783 float MCq22ev = (MCtpcQQ[2]-M)/W_2;
1784 float MCq22Pev = (MCtpcQQp[2]-Mp)/Wp_2;
1785 float MCq22Nev = (MCtpcQQn[2]-Mn)/Wn_2;
1786 float MCq22gap0ev = (MCtpcXl2*MCtpcXr2+MCtpcYl2*MCtpcYr2)/(MCtpcXl0*MCtpcXr0);
1787 float MCq22gap0Pev = (MCtpcXpl[2]*MCtpcXpr[2]+MCtpcYpl[2]*MCtpcYpr[2])/(MCtpcXpl[0]*MCtpcXpr[0]);
1788 float MCq22gap0Nev = (MCtpcXnl[2]*MCtpcXnr[2]+MCtpcYnl[2]*MCtpcYnr[2])/(MCtpcXnl[0]*MCtpcXnr[0]);
1789 float MCq22gap1ev = (MCtpcXl2o*MCtpcXr2o+MCtpcYl2o*MCtpcYr2o)/(MCtpcXl0o*MCtpcXr0o);
1790 float MCq22gap1Pev = (MCtpcXplo[2]*MCtpcXpro[2]+MCtpcYplo[2]*MCtpcYpro[2])/(MCtpcXplo[0]*MCtpcXpro[0]);
1791 float MCq22gap1Nev = (MCtpcXnlo[2]*MCtpcXnro[2]+MCtpcYnlo[2]*MCtpcYnro[2])/(MCtpcXnlo[0]*MCtpcXnro[0]);
1792 h_MCq22_cent->Fill(cent,MCq22ev);
1793 h_MCq22_centP->Fill(cent,MCq22Pev);
1794 h_MCq22_centN->Fill(cent,MCq22Nev);
1795 h_MCAq22_cent->Fill(cent,MCq22ev*MCqasymm);
1796 h_MCAq22_centP->Fill(cent,MCq22Pev*MCqasymm);
1797 h_MCAq22_centN->Fill(cent,MCq22Nev*MCqasymm);
1798 h_MCq22gap0_cent->Fill(cent,MCq22gap0ev);
1799 h_MCq22gap0_centP->Fill(cent,MCq22gap0Pev);
1800 h_MCq22gap0_centN->Fill(cent,MCq22gap0Nev);
1801 h_MCq22gap1_cent->Fill(cent,MCq22gap1ev);
1802 h_MCq22gap1_centP->Fill(cent,MCq22gap1Pev);
1803 h_MCq22gap1_centN->Fill(cent,MCq22gap1Nev);
1804
1805
1806
1807 if(debug>0) cout<<"there are "<<ntrkMC<<" Monte Carlo tracks within the cuts in this event"<<endl;
1808
1809 // --- fsecond MC track loop
1810 for(int itrkMC = 0; itrkMC<d_ntrkMC; itrkMC++)
1811 {
1812 //AliVParticle *trackMC = fMC->GetTrack(itrkMC);
1813 AliAODMCParticle *trackMC = (AliAODMCParticle *)fMC->GetTrack(itrkMC);
1814 if(!trackMC)
1815 {
1816 if(debug>0) cout<<"ERROR: Could not receive track "<<itrkMC<<" (mc loop)"<<endl;
1817 continue;
1818 }
1819
1820 float pt = trackMC->Pt();
1821
1822 float phi = trackMC->Phi();
1823 float eta = trackMC->Eta();
1824 int charge = trackMC->Charge();
1825 if(charge==0) continue; // new
1826 charge /= 3; // new
1827 bool pos = charge>0.0;
1828 bool neg = charge<0.0;
1829
1830
1831 // apply kinematic cuts
1832 if(fabs(eta)>outeta) continue;
1833 if(pt<ptmin||pt>ptmax) continue;
1834
1835 //if(charge!=-3&&charge!=+3) continue;// x3 by convention
1836 if(!trackMC->IsPrimary()) continue;
1837 if(!trackMC->IsPhysicalPrimary()) continue;
1838
1839 //if(abs(p0->GetPdgCode())==11) continue; //electrons
1840 //if(abs(p0->GetPdgCode())==13) continue; //electrons
1841
1842
1843 float MCtrkXpl[hmax], MCtrkYpl[hmax], MCtrkXpr[hmax], MCtrkYpr[hmax];
1844 float MCtrkXnl[hmax], MCtrkYnl[hmax], MCtrkXnr[hmax], MCtrkYnr[hmax];
1845 for(int i=0; i<hmax;i++)
1846 {
1847 MCtrkXpl[i] = 0.0;
1848 MCtrkYpl[i] = 0.0;
1849 MCtrkXpr[i] = 0.0;
1850 MCtrkYpr[i] = 0.0;
1851 MCtrkXnl[i] = 0.0;
1852 MCtrkYnl[i] = 0.0;
1853 MCtrkXnr[i] = 0.0;
1854 MCtrkYnr[i] = 0.0;
1855 }
1856 // very similar to above, except += changed to =
1857 if(pt>ptmin&&pt<ptmax)
1858 {
1859 if(eta>-outeta&&eta<-excleta)
1860 {
1861 // assign values to array elements 1..8
1862 for(int i=1; i<hmax;i++)
1863 {
1864 if(charge>0)
1865 {
1866 MCtrkXpl[i] = cos(i*phi);
1867 MCtrkYpl[i] = sin(i*phi);
1868 }
1869 if(charge<0)
1870 {
1871 MCtrkXnl[i] = cos(i*phi);
1872 MCtrkYnl[i] = sin(i*phi);
1873 }
1874 }
1875 // assign values to array element 0
1876 if(charge>0)
1877 {
1878 MCtrkXpl[0] = 1.0;
1879 MCtrkYpl[0] = pt;
1880 }
1881 if(charge<0)
1882 {
1883 MCtrkXnl[0] = 1.0;
1884 MCtrkYnl[0] = pt;
1885 }
1886 }
1887 //----------------------------
1888 if(eta>excleta&&eta<outeta) // right half of TPC
1889 {
1890 for(int i=1; i<hmax; i++)
1891 {
1892 if(charge>0)
1893 {
1894 MCtrkXpr[i] = cos(i*phi);
1895 MCtrkYpr[i] = sin(i*phi);
1896 }
1897 if(charge<0)
1898 {
1899 MCtrkXnr[i] = cos(i*phi);
1900 MCtrkYnr[i] = sin(i*phi);
1901 }
1902 }
1903 if(charge>0)
1904 {
1905 MCtrkXpr[0] = 1.0;
1906 MCtrkYpr[0] = pt;
1907 }
1908 if(charge<0)
1909 {
1910 MCtrkXnr[0] = 1.0;
1911 MCtrkYnr[0] = pt;
1912 }
1913 } // end right half of TPC
1914 } // end pT selection
1915
1916 // COME BACK HERE
1917 float MCtrkX[9], MCtrkY[9];
1918 float MCtrkXp[9], MCtrkYp[9]; // pos
1919 float MCtrkXn[9], MCtrkYn[9]; // neg
1920 for(int i=0; i<6; i++)
1921 {
1922 MCtrkX[i]=MCtrkXpl[i]+MCtrkXnl[i]+MCtrkXpr[i]+MCtrkXnr[i];
1923 MCtrkY[i]=MCtrkYpl[i]+MCtrkYnl[i]+MCtrkYpr[i]+MCtrkYnr[i];
1924 // pos
1925 MCtrkXp[i]=MCtrkXpl[i]+MCtrkXpr[i];
1926 MCtrkYp[i]=MCtrkYpl[i]+MCtrkYpr[i];
1927 // neg
1928 MCtrkXn[i]=MCtrkXnl[i]+MCtrkXnr[i];
1929 MCtrkYn[i]=MCtrkYnl[i]+MCtrkYnr[i];
1930 }
1931
1932 // sanity checks to prevent division by zero and other problems
1933 if(MCtrkX[0]<1) continue;
1934 if(MCtpcX[0]<2) continue;
1935 if(MCtpcXp[0]<2) continue;
1936 if(MCtpcXn[0]<2) continue;
1937
1938 // --- calculate differential cumulants from Q-vector components
1939 // --- reference flow is from whole TPC
1940 float MCdiffq22ev = ((MCtrkX[2]*MCtpcX[2])+(MCtrkY[2]*MCtpcY[2])-1)/(MCtpcX[0]-1);
1941 float MCdiffq22Pev = ((MCtrkXp[2]*MCtpcX[2])+(MCtrkYp[2]*MCtpcY[2])-1)/(MCtpcX[0]-1);
1942 float MCdiffq22Nev = ((MCtrkXn[2]*MCtpcX[2])+(MCtrkYn[2]*MCtpcY[2])-1)/(MCtpcX[0]-1);
1943 // --- third harmonic
1944 float MCdiffq32ev = ((MCtrkX[3]*MCtpcX[3])+(MCtrkY[3]*MCtpcY[3])-1)/(MCtpcX[0]-1);
1945 float MCdiffq32Pev = ((MCtrkXp[3]*MCtpcX[3])+(MCtrkYp[3]*MCtpcY[3])-1)/(MCtpcX[0]-1);
1946 float MCdiffq32Nev = ((MCtrkXn[3]*MCtpcX[3])+(MCtrkYn[3]*MCtpcY[3])-1)/(MCtpcX[0]-1);
1947 // --- fourth harmonic
1948 float MCdiffq42ev = ((MCtrkX[4]*MCtpcX[4])+(MCtrkY[4]*MCtpcY[4])-1)/(MCtpcX[0]-1);
1949 float MCdiffq42Pev = ((MCtrkXp[4]*MCtpcX[4])+(MCtrkYp[4]*MCtpcY[4])-1)/(MCtpcX[0]-1);
1950 float MCdiffq42Nev = ((MCtrkXn[4]*MCtpcX[4])+(MCtrkYn[4]*MCtpcY[4])-1)/(MCtpcX[0]-1);
1951
1952
1953 if(donested)
1954 {
1955 for(int i=0; i<d_ntrkMC; i++)
1956 {
1957 // ----------------------------------------------------------------------------------------
1958 // --- want delta eta and delta pt dependence of differential cumulant weighted with charge
1959 // --- need nested track loop to get eta1 and eta3
1960 // ----------------------------------------------------------------------------------------
1961 // make sure eta is in range, -99 should be autorejected by this
1962 //if(MCeta3[i]==-99) continue;
1963 if(MCeta3[i]<-outeta||MCeta3[i]>outeta)
1964 {
1965 //cout<<"MCeta3 out of range!!! "<<MCeta3[i]<<endl;
1966 continue;
1967 }
1968 if(eta<-outeta||eta>outeta)
1969 {
1970 cout<<"eta out of range!!! "<<eta<<endl;
1971 continue;
1972 }
1973 if(MCpt3[i]<ptmin||MCpt3[i]>ptmax)
1974 {
1975 //cout<<"MCpt3 out of range!!! "<<MCpt3[i]<<endl;
1976 continue;
1977 }
1978 if(pt<ptmin||pt>ptmax)
1979 {
1980 cout<<"pt out of range!!!"<<endl;
1981 continue;
1982 }
1983 // charge==0 should be rejected by the above cut???
1984 if(MCcharge3[i]!=1&&MCcharge3[i]!=-1)
1985 {
1986 //cout<<"WTF!"<<endl;
1987 continue;
1988 }
1989 // VERY IMPORTANT remove auto-correlations with 3rd particle
1990 if(i==itrkMC) continue;
1991
1992 float DETA = eta-MCeta3[i];
1993 //float DPHI = phi-MCphi3[i];
1994 //float DPT = pt-MCpt3[i];
1995
1996 // ---
1997
1998 hMC_AT_X_deta->Fill(DETA,MCcharge3[i]);
1999 if(pos) hMC_AT_X_detaP->Fill(DETA,MCcharge3[i]);
2000 if(neg) hMC_AT_X_detaN->Fill(DETA,MCcharge3[i]);
2001 if(pos) hMC_diffq22_X_detaP->Fill(DETA,MCdiffq22Pev);
2002 if(neg) hMC_diffq22_X_detaN->Fill(DETA,MCdiffq22Nev);
2003 if(pos) hMC_diffq32_X_detaP->Fill(DETA,MCdiffq32Pev);
2004 if(neg) hMC_diffq32_X_detaN->Fill(DETA,MCdiffq32Nev);
2005 if(pos) hMC_diffq42_X_detaP->Fill(DETA,MCdiffq42Pev);
2006 if(neg) hMC_diffq42_X_detaN->Fill(DETA,MCdiffq42Nev);
2007 if(pos) hMC_ATdiffq22_X_detaP->Fill(DETA,MCdiffq22Pev*MCcharge3[i]);
2008 if(neg) hMC_ATdiffq22_X_detaN->Fill(DETA,MCdiffq22Nev*MCcharge3[i]);
2009 if(pos) hMC_ATdiffq32_X_detaP->Fill(DETA,MCdiffq32Pev*MCcharge3[i]);
2010 if(neg) hMC_ATdiffq32_X_detaN->Fill(DETA,MCdiffq32Nev*MCcharge3[i]);
2011 if(pos) hMC_ATdiffq42_X_detaP->Fill(DETA,MCdiffq42Pev*MCcharge3[i]);
2012 if(neg) hMC_ATdiffq42_X_detaN->Fill(DETA,MCdiffq42Nev*MCcharge3[i]);
2013
2014 } // nested track loop
2015
2016 } // check on donested
2017
2018 } // end of second MC track loop
2019
2020 } // check on existence of MC
2021
2022
2023
2024 // ------------------------------- //
2025 // --- Now looping over tracks --- //
2026 // ------------------------------- //
2027
2028 // --- track loop for mapping matrix
2029 TExMap *trackMap = new TExMap();
2030 for(int itrk=0; itrk<d_ntrk; itrk++)
2031 {
2032 AliAODTrack *track = fAOD->GetTrack(itrk);
2033 if(!track)
2034 {
2035 if(debug>0) cout<<"ERROR: Could not retrieve AODtrack "<<itrk<<endl;
2036 continue;
2037 }
2038 int gid = track->GetID();
2039 if(track->TestFilterBit(fbit)) trackMap->Add(gid,itrk);
2040 }
2041
2042
2043
2044 // variables for 3rd particle
2045 float pt3[d_ntrk];
2046 float eta3[d_ntrk];
2047 float phi3[d_ntrk];//new...
2048 int charge3[d_ntrk];
2049 // initial variables to out of range values
2050 // for cases when the track loop skips certain tracks
2051 for(int i=0; i<d_ntrk; i++)
2052 {
2053 pt3[i] = -99;
2054 eta3[i] = -99;
2055 phi3[i] = -99;
2056 charge3[i] = 0;
2057 }
2058
2059 int badbit = 0; // counter for number of tracks with wrong filter bit
2060 // --- main track loop
2061 for(int itrk = 0; itrk<d_ntrk; itrk++)
2062 {
2063 AliAODTrack *track = fAOD->GetTrack(itrk);
2064 if(!track)
2065 {
2066 if(debug>0) cout<<"ERROR: Could not retrieve AODtrack "<<itrk<<endl;
2067 continue;
2068 }
2069
2070 if(!track->TestFilterBit(fbit))
2071 {
2072 badbit++; // count tracks with wrong filter bit
2073 if(debug>15) cout<<"wrong filter bit for track "<<itrk<<endl;
2074 continue;
2075 }
2076
2077 int gid = track->GetID();
2078 AliAODTrack *PIDtrack;
2079 if(gid>=0) PIDtrack = track;
2080 else PIDtrack = fAOD->GetTrack(trackMap->GetValue(-1-gid));
2081
2082
2083 // if(debug>15)
2084 // {
2085 // cout<<"track index is "<<itrk<<endl;
2086 // cout<<"track ID is "<<gid<<endl;
2087 // }
2088
2089 // --- get track variables
2090 float pt = track->Pt();
2091 float phi = track->Phi();
2092 float eta = track->Eta();
2093 int charge = track->Charge();
2094 int ncls = track->GetTPCNcls();
2095 float dedx = track->GetTPCsignal();
2096
2097 float Tdcaxy = track->DCA();
2098 float Tdcaz = track->ZAtDCA();
2099
2100
2101 float dcaxy = -999;
2102 float dcax = -999;
2103 float dcay = -999;
2104 float dcaz = -999;
2105
2106 double r[3];
2107 bool dcaflag = track->GetXYZ(r);
2108 //cout<<"fbit is "<<fbit<<" GetXYZ is "<<dcaflag<<endl;
2109 //printf("fbit is %d GetXYZ is %d PropagateToDCA is %d Tdcaz is %f pt is %f \n\n",fbit,dcaflag,proptodca,Tdcaz,pt);
2110
2111 double DCA[2]; // dca
2112 double COV[3]; // covariance
2113 bool proptodca = track->PropagateToDCA(fVtx,mag,100.0,DCA,COV);
2114 if(!proptodca&&debug>0) cout<<"No DCACOV for you!"<<endl;
2115
2116 if(dcaflag)
2117 {
2118 dcaxy = r[0];
2119 dcaz = r[1];
2120 if(debug>5) cout<<"GetXYZ is true, filter bit is "<<fbit<<endl;
2121 }
2122 else
2123 {
2124 dcax = r[0] - eventX;
2125 dcay = r[1] - eventY;
2126 dcaz = r[2] - eventZ;
2127 // --- need the sign convention for dcaxy...
2128 dcaxy = sqrt(dcax*dcax+dcay*dcay);
2129 //if((float)dcaxy!=(float)fabs((float)DCA[0])&&debug>4) cout<<"hmm... "<<dcaxy<<" "<<DCA[0]<<endl;
2130 // --- set dcaxy to value from PropagateToDCA to get correct sign
2131 dcaxy = DCA[0];
2132 // --- dcaz on the other hand is unambiguous
2133 //if(dcaz!=(float)DCA[1]&&debug>4) cout<<"hmm... "<<dcaz<<" "<<DCA[1]<<endl;
2134 if(debug>5) cout<<"GetXYZ is false, filter bit is "<<fbit<<endl;
2135 }
2136
2137 if(debug>5)
2138 {
2139 cout<<"r[0] is "<<r[0]<<" ";
2140 cout<<"r[1] is "<<r[1]<<" ";
2141 cout<<"r[2] is "<<r[2]<<endl;
2142 cout<<"eventX is "<<eventX<<" ";
2143 cout<<"eventY is "<<eventY<<" ";
2144 cout<<"eventZ is "<<eventZ<<endl;
2145 cout<<"Tdcaxy is "<<Tdcaxy<<" and dcaxy is "<<dcaxy<<" and DCACOV xy is "<<DCA[0]<<endl;
2146 cout<<"Tdcaz is "<<Tdcaz<<" and dcaz is "<<dcaz<<" and DCACOV z is "<<DCA[1]<<endl;
2147 }
2148
2149 if((fbit==128||fbit==272)&&!dcaflag)
2150 {
2151 dcaxy = Tdcaxy;
2152 dcaz = Tdcaz;
2153 if(debug>4)
2154 {
2155 cout<<"GetXYZ false but should be true, flipping values"<<endl;
2156 cout<<"Tdcaxy is "<<Tdcaxy<<" and dcaxy is "<<dcaxy<<" and DCACOV xy is "<<DCA[0]<<endl;
2157 cout<<"Tdcaz is "<<Tdcaz<<" and dcaz is "<<dcaz<<" and DCACOV z is "<<DCA[1]<<endl;
2158 }
2159 }
2160
2161 // --- some diagnostic histograms for tracks
2162 fHistPt->Fill(pt);
2163 fHistPhi->Fill(phi);
2164 fHistEta->Fill(eta);
2165 fHistCharge->Fill(charge);
2166 fHistTPCncls->Fill(ncls);
2167 fHistDedx->Fill(dedx);
2168 fHistDCAxy->Fill(dcaxy);
2169 fHistDCAz->Fill(dcaz);
2170
2171 if(charge>0)
2172 {
2173 fHistPosPt->Fill(pt);
2174 fHistPosPhi->Fill(phi);
2175 fHistPosEta->Fill(eta);
2176 }
2177 if(charge<0)
2178 {
2179 fHistNegPt->Fill(pt);
2180 fHistNegPhi->Fill(phi);
2181 fHistNegEta->Fill(eta);
2182 }
2183 fProfMeanChargePt->Fill(pt,charge);
2184 fProfMeanChargeEta->Fill(eta,charge);
2185
2186 // --- some extra verbose stuff
2187 if(debug>8&&(pt<ptmin||pt>ptmax)) cout<<"pt = "<<pt<<" out of range for Q-vectors"<<endl;
2188 if(debug>9) cout<<"number of tpc clusters is "<<ncls<<endl;
2189
2190 // --- track cut on number of TPC clusters
2191 if(ncls<nclscut)
2192 {
2193 if(debug>5) cout<<"number of tpc clusters is too low "<<ncls<<endl;
2194 continue;
2195 }
2196
2197 // --- track cut on dca
2198 if(fabs(dcaxy)>dcacutxy||fabs(dcaz)>dcacutz)
2199 {
2200 if(debug>7)
2201 {
2202 cout<<"dca out of range..."<<endl;
2203 cout<<"dcaxy cut is "<<dcacutxy<<endl;
2204 cout<<"dcaz cut is "<<dcacutz<<endl;
2205 cout<<"dcaxy is "<<dcaxy<<endl;
2206 cout<<"dcaz is "<<dcaz<<endl;
2207 }
2208 if(dodcacuts) continue; // problems depending on filter bit...
2209 }
2210 // Prabhat
2211 // fHistPhi->Fill(phi);
2212 // fHistEta->Fill(eta);
2213 fHistDCAxyAfter->Fill(dcaxy);
2214 fHistDCAzAfter->Fill(dcaz);
2215
2d95bc0f 2216 float nsigmapion = 0.;
2217 float nsigmakaon = 0.;
2218 float nsigmaprot = 0.;
2219 float nsigmaelec = 0.;
2220
2221 if(!doMC){
2222 nsigmapion = fPID->NumberOfSigmasTPC(PIDtrack,(AliPID::EParticleType)AliPID::kPion);
2223 nsigmakaon = fPID->NumberOfSigmasTPC(PIDtrack,(AliPID::EParticleType)AliPID::kKaon);
2224 nsigmaprot = fPID->NumberOfSigmasTPC(PIDtrack,(AliPID::EParticleType)AliPID::kProton);
2225 nsigmaelec = fPID->NumberOfSigmasTPC(PIDtrack,(AliPID::EParticleType)AliPID::kElectron);
2226 }
545d70ee 2227
2228 bool isPion = fabs(nsigmapion) <= nspid;
2229 bool isKaon = fabs(nsigmakaon) <= nspid;
2230 bool isProt = fabs(nsigmaprot) <= nspid;
2231 bool isElec = fabs(nsigmaelec) <= nspid;
2232
2233 bool isPionH = isPion && (!isKaon) && (!isProt);
2234 bool isProtH = (!isPion) && (!isKaon) && isProt;
2235 bool isPionL = isPion && (!isKaon) && (!isProt) && (!isElec);
2236 bool isProtL = (!isPion) && (!isKaon) && isProt && (!isElec);
2237
2238 if(isPion) fHistPtPion->Fill(pt);
2239 if(isPionH) fHistPtPionH->Fill(pt);
2240 if(isPionL) fHistPtPionL->Fill(pt);
2241 fHistNsigmaPion->Fill(nsigmapion);
2242
2243 if(isProt) fHistPtProt->Fill(pt);
2244 if(isProtH) fHistPtProtH->Fill(pt);
2245 if(isProtL) fHistPtProtL->Fill(pt);
2246 fHistNsigmaProt->Fill(nsigmapion);
2247
2248 // --- count track numbers
2249 ntrk++;
2250 if(charge>0) ntrkpos++;
2251 if(charge<0) ntrkneg++;
2252 if(eta<0)
2253 {
2254 ntrkL++;
2255 if(charge>0) ntrkposL++;
2256 if(charge<0) ntrknegL++;
2257 }
2258 if(eta>0)
2259 {
2260 ntrkR++;
2261 if(charge>0) ntrkposR++;
2262 if(charge<0) ntrknegR++;
2263 }
2264 h_AT_eta->Fill(eta,charge);
2265 if(mag<-4.0) h_AT_etaF1->Fill(eta,charge);
2266 if(mag>4.0) h_AT_etaF3->Fill(eta,charge);
2267 h2_AT_eta->Fill(eta,charge);
2268 if(mag<-4.0) h2_AT_etaF1->Fill(eta,charge);
2269 if(mag>4.0) h2_AT_etaF3->Fill(eta,charge);
2270 //
2271 if(mag<-4.0&&charge>0) h_eta_pos_F1->Fill(eta);
2272 if(mag<-4.0&&charge<0) h_eta_neg_F1->Fill(eta);
2273 if(mag>4.0&&charge>0) h_eta_pos_F3->Fill(eta);
2274 if(mag>4.0&&charge<0) h_eta_neg_F3->Fill(eta);
2275
2276 if(doeffcorr)
2277 {
2278 if(gRandom) delete gRandom;
2279 gRandom = new TRandom3(0);
2280 gRandom->SetSeed(0);
2281 double rand = gRandom->Rndm();
2282 int ptbin = int(pt*10);
2283 if(ptbin>49) ptbin = 49;
2284 float eff = effTPC[ptbin];
2285 if(fbit==272) eff = effHYB[ptbin];
2286 if(rand>eff)
2287 {
2288 if(debug>5) cout<<"doing random throw for efficiency correction: pt is "<<pt<<" eff is "<<eff<<" rand is "<<rand<<endl;
2289 continue;
2290 }
2291 }
2292
2293 if(pt>ptmin&&pt<ptmax&&fabs(eta)<outeta&&fabs(eta)>excleta)
2294 {
2295 cutntrk++;
2296 if(charge>0) cutntrkpos++;
2297 if(charge<0) cutntrkneg++;
2298 if(eta<0)
2299 {
2300 cutntrkL++;
2301 if(charge>0) cutntrkposL++;
2302 if(charge<0) cutntrknegL++;
2303 }
2304 if(eta>0)
2305 {
2306 cutntrkR++;
2307 if(charge>0) cutntrkposR++;
2308 if(charge<0) cutntrknegR++;
2309 }
2310 h_AT_cut_eta->Fill(eta,charge);
2311 if(mag<-4.0) h_AT_cut_etaF1->Fill(eta,charge);
2312 if(mag>4.0) h_AT_cut_etaF3->Fill(eta,charge);
2313 h2_AT_cut_eta->Fill(eta,charge);
2314 if(mag<-4.0) h2_AT_cut_etaF1->Fill(eta,charge);
2315 if(mag>4.0) h2_AT_cut_etaF3->Fill(eta,charge);
2316 //
2317 if(cent>=centlo&&cent<centhi)
2318 {
2319 if(mag<-4.0&&charge>0) h_cut_eta_pos_F1->Fill(eta);
2320 if(mag<-4.0&&charge<0) h_cut_eta_neg_F1->Fill(eta);
2321 if(mag>4.0&&charge>0) h_cut_eta_pos_F3->Fill(eta);
2322 if(mag>4.0&&charge<0) h_cut_eta_neg_F3->Fill(eta);
2323 }
2324 }
2325
2326 // ------------------------------------ //
2327 // --- now fill Q-vector components --- //
2328 // ------------------------------------ //
2329
2330 if(pt>ptmin&&pt<ptmax)
2331 {
2332 // --- outer
2333 if(eta>-outeta&&eta<-ineta)
2334 {
2335 for(int i=1; i<hmax;i++)
2336 {
2337 if(charge>0)
2338 {
2339 tpcXplo[i] += cos(i*phi);
2340 tpcYplo[i] += sin(i*phi);
2341 tpcXpl[i] += cos(i*phi);
2342 tpcYpl[i] += sin(i*phi);
2343 }
2344 if(charge<0)
2345 {
2346 tpcXnlo[i] += cos(i*phi);
2347 tpcYnlo[i] += sin(i*phi);
2348 tpcXnl[i] += cos(i*phi);
2349 tpcYnl[i] += sin(i*phi);
2350 }
2351 }
2352 if(charge>0)
2353 {
2354 tpcXplo[0] += 1.0;
2355 tpcYplo[0] += pt;
2356 tpcXpl[0] += 1.0;
2357 tpcYpl[0] += pt;
2358 }
2359 if(charge<0)
2360 {
2361 tpcXnlo[0] += 1.0;
2362 tpcYnlo[0] += pt;
2363 tpcXnl[0] += 1.0;
2364 tpcYnl[0] += pt;
2365 }
2366 } // end left half of TPC
2367 //----------------------------
2368 if(eta>ineta&&eta<outeta)
2369 {
2370 for(int i=1; i<hmax; i++)
2371 {
2372 if(charge>0)
2373 {
2374 tpcXpro[i] += cos(i*phi);
2375 tpcYpro[i] += sin(i*phi);
2376 tpcXpr[i] += cos(i*phi);
2377 tpcYpr[i] += sin(i*phi);
2378 }
2379 if(charge<0)
2380 {
2381 tpcXnro[i] += cos(i*phi);
2382 tpcYnro[i] += sin(i*phi);
2383 tpcXnr[i] += cos(i*phi);
2384 tpcYnr[i] += sin(i*phi);
2385 }
2386 }
2387 if(charge>0)
2388 {
2389 tpcXpro[0] += 1.0;
2390 tpcYpro[0] += pt;
2391 tpcXpr[0] += 1.0;
2392 tpcYpr[0] += pt;
2393 }
2394 if(charge<0)
2395 {
2396 tpcXnro[0] += 1.0;
2397 tpcYnro[0] += pt;
2398 tpcXnr[0] += 1.0;
2399 tpcYnr[0] += pt;
2400 }
2401 } // end right half of TPC
2402 // --- inner
2403 if(eta>-ineta&&eta<-excleta)
2404 {
2405 for(int i=1; i<hmax;i++)
2406 {
2407 if(charge>0)
2408 {
2409 tpcXpli[i] += cos(i*phi);
2410 tpcYpli[i] += sin(i*phi);
2411 tpcXpl[i] += cos(i*phi);
2412 tpcYpl[i] += sin(i*phi);
2413 }
2414 if(charge<0)
2415 {
2416 tpcXnli[i] += cos(i*phi);
2417 tpcYnli[i] += sin(i*phi);
2418 tpcXnl[i] += cos(i*phi);
2419 tpcYnl[i] += sin(i*phi);
2420 }
2421 }
2422 if(charge>0)
2423 {
2424 tpcXpli[0] += 1.0;
2425 tpcYpli[0] += pt;
2426 tpcXpl[0] += 1.0;
2427 tpcYpl[0] += pt;
2428 }
2429 if(charge<0)
2430 {
2431 tpcXnli[0] += 1.0;
2432 tpcYnli[0] += pt;
2433 tpcXnl[0] += 1.0;
2434 tpcYnl[0] += pt;
2435 }
2436 } // end left half of TPC
2437 //----------------------------
2438 if(eta>excleta&&eta<ineta)
2439 {
2440 for(int i=1; i<hmax; i++)
2441 {
2442 if(charge>0)
2443 {
2444 tpcXpri[i] += cos(i*phi);
2445 tpcYpri[i] += sin(i*phi);
2446 tpcXpr[i] += cos(i*phi);
2447 tpcYpr[i] += sin(i*phi);
2448 }
2449 if(charge<0)
2450 {
2451 tpcXnri[i] += cos(i*phi);
2452 tpcYnri[i] += sin(i*phi);
2453 tpcXnr[i] += cos(i*phi);
2454 tpcYnr[i] += sin(i*phi);
2455 }
2456 }
2457 if(charge>0)
2458 {
2459 tpcXpri[0] += 1.0;
2460 tpcYpri[0] += pt;
2461 tpcXpr[0] += 1.0;
2462 tpcYpr[0] += pt;
2463 }
2464 if(charge<0)
2465 {
2466 tpcXnri[0] += 1.0;
2467 tpcYnri[0] += pt;
2468 tpcXnr[0] += 1.0;
2469 tpcYnr[0] += pt;
2470 }
2471 } // end right half of TPC
2472 } // end pt selection
2473
2474 // set values for 3rd particle variables
2475 pt3[itrk] = pt;
2476 eta3[itrk] = eta;
2477 phi3[itrk] = phi;
2478 charge3[itrk] = (int)charge;
2479
2480
2481 // ------------------------------------ //
2482 // --- that's it for the track loop --- //
2483 // ------------------------------------ //
2484
2485 } // end of first track loop
2486
2487 if(debug>0)
2488 {
2489 cout<<"filter bit is "<<fbit<<endl;
2490 cout<<"number of tracks with wrong filter bit "<<badbit<<endl;
2491 cout<<"difference = ntrk - nbad = "<<d_ntrk-badbit<<endl;
2492 }
2493
2494 float tpcX[9], tpcY[9], tpcQQ[9];//, qq[9];
2495 float tpcXp[9], tpcYp[9], tpcQQp[9];//, qqp[9]; // pos
2496 float tpcXn[9], tpcYn[9], tpcQQn[9];//, qqn[9]; // neg
2497
2498 for(int i=0; i<6; i++)
2499 {
2500 tpcX[i]=tpcXpl[i]+tpcXnl[i]+tpcXpr[i]+tpcXnr[i];
2501 tpcY[i]=tpcYpl[i]+tpcYnl[i]+tpcYpr[i]+tpcYnr[i];
2502 tpcQQ[i]=tpcX[i]*tpcX[i]+tpcY[i]*tpcY[i];
2503 //qq[i]=sqrt(tpcQQ[i]/tpcX[0]);
2504 // pos
2505 tpcXp[i]=tpcXpl[i]+tpcXpr[i];
2506 tpcYp[i]=tpcYpl[i]+tpcYpr[i];
2507 tpcQQp[i]=tpcXp[i]*tpcXp[i]+tpcYp[i]*tpcYp[i];
2508 //qqp[i]=sqrt(tpcQQp[i]/tpcXp[0]);
2509 // neg
2510 tpcXn[i]=tpcXnl[i]+tpcXnr[i];
2511 tpcYn[i]=tpcYnl[i]+tpcYnr[i];
2512 tpcQQn[i]=tpcXn[i]*tpcXn[i]+tpcYn[i]*tpcYn[i];
2513 //qqn[i]=sqrt(tpcQQn[i]/tpcXn[0]);
2514 }
2515
2516
2517
2518
2519
2520 float M = tpcX[0];
2521 float W_2 = M*(M-1);
2522 float Mp = tpcXp[0];
2523 float Wp_2 = Mp*(Mp-1);
2524 float Mn = tpcXn[0];
2525 float Wn_2 = Mn*(Mn-1);
2526
2527 // prevent division by zero
2528 if(M<2) return;
2529 if(Mp<2) return;
2530 if(Mn<2) return;
2531
2532 if(ntrkL<1) return;
2533 if(ntrkR<1) return;
2534
2535 float qasymm = (float)(ntrkpos-ntrkneg)/ntrk;
2536 float qasymmL = (float)(ntrkposL-ntrknegL)/ntrkL;
2537 float qasymmR = (float)(ntrkposR-ntrknegR)/ntrkR;
2538
2539 if(doacuts)
2540 {
2541 qasymm = (float)(cutntrkpos-cutntrkneg)/cutntrk;
2542 qasymmL = (float)(cutntrkposL-cutntrknegL)/cutntrkL;
2543 qasymmR = (float)(cutntrkposR-cutntrknegR)/cutntrkR;
2544 }
2545
2546 h_A_cent->Fill(cent,qasymm);
2547 if(mag<-4.0) h_A_centF1->Fill(cent,qasymm);
2548 if(mag>4.0) h_A_centF3->Fill(cent,qasymm);
2549 h2_A_cent->Fill(cent,qasymm);
2550 if(mag<-4.0) h2_A_centF1->Fill(cent,qasymm);
2551 if(mag>4.0) h2_A_centF3->Fill(cent,qasymm);
2552
2553 float Xl2 = tpcXnl[2]+tpcXpl[2];
2554 float Yl2 = tpcYnl[2]+tpcYpl[2];
2555 float Xr2 = tpcXnr[2]+tpcXpr[2];
2556 float Yr2 = tpcYnr[2]+tpcYpr[2];
2557 float Xl3 = tpcXnl[3]+tpcXpl[3];
2558 float Yl3 = tpcYnl[3]+tpcYpl[3];
2559 float Xr3 = tpcXnr[3]+tpcXpr[3];
2560 float Yr3 = tpcYnr[3]+tpcYpr[3];
2561 float Xl4 = tpcXnl[4]+tpcXpl[4];
2562 float Yl4 = tpcYnl[4]+tpcYpl[4];
2563 float Xr4 = tpcXnr[4]+tpcXpr[4];
2564 float Yr4 = tpcYnr[4]+tpcYpr[4];
2565 float Xl0 = tpcXnl[0]+tpcXpl[0];
2566 float Xr0 = tpcXnr[0]+tpcXpr[0];
2567
2568 float Xl2o = tpcXnlo[2]+tpcXplo[2];
2569 float Yl2o = tpcYnlo[2]+tpcYplo[2];
2570 float Xr2o = tpcXnro[2]+tpcXpro[2];
2571 float Yr2o = tpcYnro[2]+tpcYpro[2];
2572 float Xl3o = tpcXnlo[3]+tpcXplo[3];
2573 float Yl3o = tpcYnlo[3]+tpcYplo[3];
2574 float Xr3o = tpcXnro[3]+tpcXpro[3];
2575 float Yr3o = tpcYnro[3]+tpcYpro[3];
2576 float Xl4o = tpcXnlo[4]+tpcXplo[4];
2577 float Yl4o = tpcYnlo[4]+tpcYplo[4];
2578 float Xr4o = tpcXnro[4]+tpcXpro[4];
2579 float Yr4o = tpcYnro[4]+tpcYpro[4];
2580 float Xl0o = tpcXnlo[0]+tpcXplo[0];
2581 float Xr0o = tpcXnro[0]+tpcXpro[0];
2582
2583 // --- q vector components for recentering if needed
2584 h_X22_cent->Fill(cent,(tpcX[2]/tpcX[0]));
2585 h_X22_centP->Fill(cent,(tpcXp[2]/tpcXp[0]));
2586 h_X22_centN->Fill(cent,(tpcXn[2]/tpcXn[0]));
2587 h_Y22_cent->Fill(cent,(tpcY[2]/tpcX[0]));
2588 h_Y22_centP->Fill(cent,(tpcYp[2]/tpcXp[0]));
2589 h_Y22_centN->Fill(cent,(tpcYn[2]/tpcXn[0]));
2590 // --- outer left
2591 h_X22_lo_cent->Fill(cent,(tpcX[2]/tpcX[0]));
2592 h_X22_lo_centP->Fill(cent,(tpcXplo[2]/tpcXplo[0]));
2593 h_X22_lo_centN->Fill(cent,(tpcXnlo[2]/tpcXnlo[0]));
2594 h_Y22_lo_cent->Fill(cent,(tpcY[2]/tpcX[0]));
2595 h_Y22_lo_centP->Fill(cent,(tpcYplo[2]/tpcXplo[0]));
2596 h_Y22_lo_centN->Fill(cent,(tpcYnlo[2]/tpcXnlo[0]));
2597 // --- inner left
2598 h_X22_li_cent->Fill(cent,(tpcX[2]/tpcX[0]));
2599 h_X22_li_centP->Fill(cent,(tpcXpli[2]/tpcXpli[0]));
2600 h_X22_li_centN->Fill(cent,(tpcXnli[2]/tpcXnli[0]));
2601 h_Y22_li_cent->Fill(cent,(tpcY[2]/tpcX[0]));
2602 h_Y22_li_centP->Fill(cent,(tpcYpli[2]/tpcXpli[0]));
2603 h_Y22_li_centN->Fill(cent,(tpcYnli[2]/tpcXnli[0]));
2604 // --- inner right
2605 h_X22_ri_cent->Fill(cent,(tpcX[2]/tpcX[0]));
2606 h_X22_ri_centP->Fill(cent,(tpcXpri[2]/tpcXpri[0]));
2607 h_X22_ri_centN->Fill(cent,(tpcXnri[2]/tpcXnri[0]));
2608 h_Y22_ri_cent->Fill(cent,(tpcY[2]/tpcX[0]));
2609 h_Y22_ri_centP->Fill(cent,(tpcYpri[2]/tpcXpri[0]));
2610 h_Y22_ri_centN->Fill(cent,(tpcYnri[2]/tpcXnri[0]));
2611 // --- outer right
2612 h_X22_ro_cent->Fill(cent,(tpcX[2]/tpcX[0]));
2613 h_X22_ro_centP->Fill(cent,(tpcXpro[2]/tpcXpro[0]));
2614 h_X22_ro_centN->Fill(cent,(tpcXnro[2]/tpcXnro[0]));
2615 h_Y22_ro_cent->Fill(cent,(tpcY[2]/tpcX[0]));
2616 h_Y22_ro_centP->Fill(cent,(tpcYpro[2]/tpcXpro[0]));
2617 h_Y22_ro_centN->Fill(cent,(tpcYnro[2]/tpcXnro[0]));
2618
2619 float q22ev = (tpcQQ[2]-M)/W_2;
2620 float q22Pev = (tpcQQp[2]-Mp)/Wp_2;
2621 float q22Nev = (tpcQQn[2]-Mn)/Wn_2;
2622 float q22gap0ev = (Xl2*Xr2+Yl2*Yr2)/(Xl0*Xr0);
2623 float q22gap0Pev = (tpcXpl[2]*tpcXpr[2]+tpcYpl[2]*tpcYpr[2])/(tpcXpl[0]*tpcXpr[0]);
2624 float q22gap0Nev = (tpcXnl[2]*tpcXnr[2]+tpcYnl[2]*tpcYnr[2])/(tpcXnl[0]*tpcXnr[0]);
2625 float q22gap1ev = (Xl2o*Xr2o+Yl2o*Yr2o)/(Xl0o*Xr0o);
2626 float q22gap1Pev = (tpcXplo[2]*tpcXpro[2]+tpcYplo[2]*tpcYpro[2])/(tpcXplo[0]*tpcXpro[0]);
2627 float q22gap1Nev = (tpcXnlo[2]*tpcXnro[2]+tpcYnlo[2]*tpcYnro[2])/(tpcXnlo[0]*tpcXnro[0]);
2628 float q22SLev = (Xl2*Xr2+Yl2*Yr2)/(Xl0*Xr0);
2629 float q22SLPev = (tpcXpl[2]*Xr2+tpcYpl[2]*Yr2)/(tpcXpl[0]*Xr0);
2630 float q22SLNev = (tpcXnl[2]*Xr2+tpcYnl[2]*Yr2)/(tpcXnl[0]*Xr0);
2631 h_q22_cent->Fill(cent,q22ev);
2632 h_q22_centP->Fill(cent,q22Pev);
2633 h_q22_centN->Fill(cent,q22Nev);
2634 h_q22gap0_cent->Fill(cent,q22gap0ev);
2635 h_q22gap0_centP->Fill(cent,q22gap0Pev);
2636 h_q22gap0_centN->Fill(cent,q22gap0Nev);
2637 h_q22gap1_cent->Fill(cent,q22gap1ev);
2638 h_q22gap1_centP->Fill(cent,q22gap1Pev);
2639 h_q22gap1_centN->Fill(cent,q22gap1Nev);
2640
2641 float q32ev = (tpcQQ[3]-M)/W_2;
2642 float q32Pev = (tpcQQp[3]-Mp)/Wp_2;
2643 float q32Nev = (tpcQQn[3]-Mn)/Wn_2;
2644 float q32gap0ev = (Xl3*Xr3+Yl3*Yr3)/(Xl0*Xr0);
2645 float q32gap0Pev = (tpcXpl[3]*tpcXpr[3]+tpcYpl[3]*tpcYpr[3])/(tpcXpl[0]*tpcXpr[0]);
2646 float q32gap0Nev = (tpcXnl[3]*tpcXnr[3]+tpcYnl[3]*tpcYnr[3])/(tpcXnl[0]*tpcXnr[0]);
2647 float q32gap1ev = (Xl3o*Xr3o+Yl3o*Yr3o)/(Xl0o*Xr0o);
2648 float q32gap1Pev = (tpcXplo[3]*tpcXpro[3]+tpcYplo[3]*tpcYpro[3])/(tpcXplo[0]*tpcXpro[0]);
2649 float q32gap1Nev = (tpcXnlo[3]*tpcXnro[3]+tpcYnlo[3]*tpcYnro[3])/(tpcXnlo[0]*tpcXnro[0]);
2650 float q32SLev = (Xl3*Xr3+Yl3*Yr3)/(Xl0*Xr0);
2651 float q32SLPev = (tpcXpl[3]*Xr3+tpcYpl[3]*Yr3)/(tpcXpl[0]*Xr0);
2652 float q32SLNev = (tpcXnl[3]*Xr3+tpcYnl[3]*Yr3)/(tpcXnl[0]*Xr0);
2653
2654 float q42ev = (tpcQQ[4]-M)/W_2;
2655 float q42Pev = (tpcQQp[4]-Mp)/Wp_2;
2656 float q42Nev = (tpcQQn[4]-Mn)/Wn_2;
2657 float q42gap0ev = (Xl4*Xr4+Yl4*Yr4)/(Xl0*Xr0);
2658 float q42gap0Pev = (tpcXpl[4]*tpcXpr[4]+tpcYpl[4]*tpcYpr[4])/(tpcXpl[0]*tpcXpr[0]);
2659 float q42gap0Nev = (tpcXnl[4]*tpcXnr[4]+tpcYnl[4]*tpcYnr[4])/(tpcXnl[0]*tpcXnr[0]);
2660 float q42gap1ev = (Xl4o*Xr4o+Yl4o*Yr4o)/(Xl0o*Xr0o);
2661 float q42gap1Pev = (tpcXplo[4]*tpcXpro[4]+tpcYplo[4]*tpcYpro[4])/(tpcXplo[0]*tpcXpro[0]);
2662 float q42gap1Nev = (tpcXnlo[4]*tpcXnro[4]+tpcYnlo[4]*tpcYnro[4])/(tpcXnlo[0]*tpcXnro[0]);
2663 float q42SLev = (Xl4*Xr4+Yl4*Yr4)/(Xl0*Xr0);
2664 float q42SLPev = (tpcXpl[4]*Xr4+tpcYpl[4]*Yr4)/(tpcXpl[0]*Xr0);
2665 float q42SLNev = (tpcXnl[4]*Xr4+tpcYnl[4]*Yr4)/(tpcXnl[0]*Xr0);
2666
2667
2668 // 3p correlator
2669 float q23ev = calc3(tpcX[1],tpcY[1],tpcX[2],tpcY[2],M);
2670 float q23Pev = calc3(tpcXp[1],tpcYp[1],tpcXp[2],tpcYp[2],Mp);
2671 float q23Nev = calc3(tpcXn[1],tpcYn[1],tpcXn[2],tpcYn[2],Mn);
2672 // 4p correlator
2673 float q24ev = calc4(tpcX[2],tpcY[2],tpcX[4],tpcY[4],M);
2674 float q24Pev = calc4(tpcXp[2],tpcYp[2],tpcXp[4],tpcYp[4],Mp);
2675 float q24Nev = calc4(tpcXn[2],tpcYn[2],tpcXn[4],tpcYn[4],Mn);
2676
2677
2678
2679 h_q22_cent->Fill(cent,q22ev);
2680 h_q22_centP->Fill(cent,q22Pev);
2681 h_q22_centN->Fill(cent,q22Nev);
2682 h_q23_cent->Fill(cent,q23ev);
2683 h_q23_centP->Fill(cent,q23Pev);
2684 h_q23_centN->Fill(cent,q23Nev);
2685 h_q24_cent->Fill(cent,q24ev);
2686 h_q24_centP->Fill(cent,q24Pev);
2687 h_q24_centN->Fill(cent,q24Nev);
2688 h_q22gap0_cent->Fill(cent,q22gap0ev);
2689 h_q22gap0_centP->Fill(cent,q22gap0Pev);
2690 h_q22gap0_centN->Fill(cent,q22gap0Nev);
2691 h_q22gap1_cent->Fill(cent,q22gap1ev);
2692 h_q22gap1_centP->Fill(cent,q22gap1Pev);
2693 h_q22gap1_centN->Fill(cent,q22gap1Nev);
2694
2695 h_q32_cent->Fill(cent,q32ev);
2696 h_q32_centP->Fill(cent,q32Pev);
2697 h_q32_centN->Fill(cent,q32Nev);
2698 h_q32gap0_cent->Fill(cent,q32gap0ev);
2699 h_q32gap0_centP->Fill(cent,q32gap0Pev);
2700 h_q32gap0_centN->Fill(cent,q32gap0Nev);
2701 h_q32gap1_cent->Fill(cent,q32gap1ev);
2702 h_q32gap1_centP->Fill(cent,q32gap1Pev);
2703 h_q32gap1_centN->Fill(cent,q32gap1Nev);
2704
2705 h_q42_cent->Fill(cent,q42ev);
2706 h_q42_centP->Fill(cent,q42Pev);
2707 h_q42_centN->Fill(cent,q42Nev);
2708 h_q42gap0_cent->Fill(cent,q42gap0ev);
2709 h_q42gap0_centP->Fill(cent,q42gap0Pev);
2710 h_q42gap0_centN->Fill(cent,q42gap0Nev);
2711 h_q42gap1_cent->Fill(cent,q42gap1ev);
2712 h_q42gap1_centP->Fill(cent,q42gap1Pev);
2713 h_q42gap1_centN->Fill(cent,q42gap1Nev);
2714
2715 if(fabs(qasymm)<10.0&&fabs(qasymmL)<10.0&&fabs(qasymmR)<10.0)
2716 {
2717 h_q22qasym_cent[icent]->Fill(qasymm,q22ev);
2718 h_q22qasymP_cent[icent]->Fill(qasymm,q22Pev);
2719 h_q22qasymN_cent[icent]->Fill(qasymm,q22Nev);
2720 h_q22qasym_gap0_cent[icent]->Fill(qasymm,q22gap0ev);
2721 h_q22qasymP_gap0_cent[icent]->Fill(qasymm,q22gap0Pev);
2722 h_q22qasymN_gap0_cent[icent]->Fill(qasymm,q22gap0Nev);
2723 h_q22qasym_gap1_cent[icent]->Fill(qasymm,q22gap1ev);
2724 h_q22qasymP_gap1_cent[icent]->Fill(qasymm,q22gap1Pev);
2725 h_q22qasymN_gap1_cent[icent]->Fill(qasymm,q22gap1Nev);
2726
2727 h_q32qasym_cent[icent]->Fill(qasymm,q32ev);
2728 h_q32qasymP_cent[icent]->Fill(qasymm,q32Pev);
2729 h_q32qasymN_cent[icent]->Fill(qasymm,q32Nev);
2730 h_q32qasym_gap0_cent[icent]->Fill(qasymm,q32gap0ev);
2731 h_q32qasymP_gap0_cent[icent]->Fill(qasymm,q32gap0Pev);
2732 h_q32qasymN_gap0_cent[icent]->Fill(qasymm,q32gap0Nev);
2733 h_q32qasym_gap1_cent[icent]->Fill(qasymm,q32gap1ev);
2734 h_q32qasymP_gap1_cent[icent]->Fill(qasymm,q32gap1Pev);
2735 h_q32qasymN_gap1_cent[icent]->Fill(qasymm,q32gap1Nev);
2736
2737 h_q42qasym_cent[icent]->Fill(qasymm,q42ev);
2738 h_q42qasymP_cent[icent]->Fill(qasymm,q42Pev);
2739 h_q42qasymN_cent[icent]->Fill(qasymm,q42Nev);
2740 h_q42qasym_gap0_cent[icent]->Fill(qasymm,q42gap0ev);
2741 h_q42qasymP_gap0_cent[icent]->Fill(qasymm,q42gap0Pev);
2742 h_q42qasymN_gap0_cent[icent]->Fill(qasymm,q42gap0Nev);
2743 h_q42qasym_gap1_cent[icent]->Fill(qasymm,q42gap1ev);
2744 h_q42qasymP_gap1_cent[icent]->Fill(qasymm,q42gap1Pev);
2745 h_q42qasymN_gap1_cent[icent]->Fill(qasymm,q42gap1Nev);
2746 }
2747
2748 h_Aq22_cent->Fill(cent,q22ev*qasymm);
2749 h_Aq22_centP->Fill(cent,q22Pev*qasymm);
2750 h_Aq22_centN->Fill(cent,q22Nev*qasymm);
2751 h_Aq22_SL_cent->Fill(cent,q22SLev*qasymm);
2752 h_Aq22_SL_centP->Fill(cent,q22SLPev*qasymm);
2753 h_Aq22_SL_centN->Fill(cent,q22SLNev*qasymm);
2754 h_Aq22_SLL_cent->Fill(cent,q22SLev*qasymmL);
2755 h_Aq22_SLL_centP->Fill(cent,q22SLPev*qasymmL);
2756 h_Aq22_SLL_centN->Fill(cent,q22SLNev*qasymmL);
2757 h_Aq22_SLR_cent->Fill(cent,q22SLev*qasymmR);
2758 h_Aq22_SLR_centP->Fill(cent,q22SLPev*qasymmR);
2759 h_Aq22_SLR_centN->Fill(cent,q22SLNev*qasymmR);
2760
2761 h_Aq32_cent->Fill(cent,q32ev*qasymm);
2762 h_Aq32_centP->Fill(cent,q32Pev*qasymm);
2763 h_Aq32_centN->Fill(cent,q32Nev*qasymm);
2764 h_Aq32_SL_cent->Fill(cent,q32SLev*qasymm);
2765 h_Aq32_SL_centP->Fill(cent,q32SLPev*qasymm);
2766 h_Aq32_SL_centN->Fill(cent,q32SLNev*qasymm);
2767 h_Aq32_SLL_cent->Fill(cent,q32SLev*qasymmL);
2768 h_Aq32_SLL_centP->Fill(cent,q32SLPev*qasymmL);
2769 h_Aq32_SLL_centN->Fill(cent,q32SLNev*qasymmL);
2770 h_Aq32_SLR_cent->Fill(cent,q32SLev*qasymmR);
2771 h_Aq32_SLR_centP->Fill(cent,q32SLPev*qasymmR);
2772 h_Aq32_SLR_centN->Fill(cent,q32SLNev*qasymmR);
2773
2774 h_Aq42_cent->Fill(cent,q42ev*qasymm);
2775 h_Aq42_centP->Fill(cent,q42Pev*qasymm);
2776 h_Aq42_centN->Fill(cent,q42Nev*qasymm);
2777 h_Aq42_SL_cent->Fill(cent,q42SLev*qasymm);
2778 h_Aq42_SL_centP->Fill(cent,q42SLPev*qasymm);
2779 h_Aq42_SL_centN->Fill(cent,q42SLNev*qasymm);
2780 h_Aq42_SLL_cent->Fill(cent,q42SLev*qasymmL);
2781 h_Aq42_SLL_centP->Fill(cent,q42SLPev*qasymmL);
2782 h_Aq42_SLL_centN->Fill(cent,q42SLNev*qasymmL);
2783 h_Aq42_SLR_cent->Fill(cent,q42SLev*qasymmR);
2784 h_Aq42_SLR_centP->Fill(cent,q42SLPev*qasymmR);
2785 h_Aq42_SLR_centN->Fill(cent,q42SLNev*qasymmR);
2786
2787
2788
2789 float a1Xp[hmax], a1Yp[hmax], a2Xp[hmax], a2Yp[hmax];
2790 float a1Xn[hmax], a1Yn[hmax], a2Xn[hmax], a2Yn[hmax];
2791
2792 for(int i=0; i<hmax;i++)
2793 {
2794 a1Xp[i] = 0.0;
2795 a1Yp[i] = 0.0;
2796 a1Xn[i] = 0.0;
2797 a1Yn[i] = 0.0;
2798 a2Xn[i] = 0.0;
2799 a2Yn[i] = 0.0;
2800 a2Xp[i] = 0.0;
2801 a2Yp[i] = 0.0;
2802 }
2803
2804 // --- now second track loop
2805
2806 for(int itrk = 0; itrk<d_ntrk; itrk++)
2807 {
2808 AliAODTrack *track = fAOD->GetTrack(itrk);
2809 if(!track)
2810 {
2811 if(debug>0) cout<<"ERROR: Could not retrieve AODtrack "<<itrk<<endl;
2812 continue;
2813 }
2814
2815 if(!track->TestFilterBit(fbit))
2816 {
2817 if(debug>15) cout<<"second track loop: wrong filter bit for track "<<itrk<<endl;
2818 continue;
2819 }
2820
2821 // int gid = track->GetID();
2822 // AliAODTrack *PIDtrack;
2823 // if(gid>=0) PIDtrack = track;
2824 // else PIDtrack = fAOD->GetTrack(trackMap->GetValue(-1-gid));
2825
2826 // --- get track variables
2827 float pt = track->Pt();
2828 float phi = track->Phi();
2829 float eta = track->Eta();
2830 int charge = track->Charge();
2831 int ncls = track->GetTPCNcls();
2832 //float dedx = track->GetTPCsignal();
2833
2834 float Tdcaxy = track->DCA();
2835 float Tdcaz = track->ZAtDCA();
2836
2837
2838 float dcaxy = -999;
2839 float dcax = -999;
2840 float dcay = -999;
2841 float dcaz = -999;
2842
2843 double r[3];
2844 bool dcaflag = track->GetXYZ(r);
2845
2846 double DCA[2]; // dca
2847 double COV[3]; // covariance
2848 bool proptodca = track->PropagateToDCA(fVtx,mag,100.0,DCA,COV);
2849 if(!proptodca&&debug>0) cout<<"No DCACOV for you!"<<endl;
2850
2851 if(dcaflag)
2852 {
2853 dcaxy = r[0];
2854 dcaz = r[1];
2855 //if(debug>6) cout<<"GetXYZ is true, filter bit is "<<fbit<<endl;
2856 }
2857 else
2858 {
2859 dcax = r[0] - eventX;
2860 dcay = r[1] - eventY;
2861 dcaz = r[2] - eventZ;
2862 // --- need the sign convention for dcaxy...
2863 dcaxy = sqrt(dcax*dcax+dcay*dcay);
2864 //if((float)dcaxy!=(float)fabs((float)DCA[0])&&debug>6) cout<<"hmm... "<<dcaxy<<" "<<DCA[0]<<endl;
2865 // --- set dcaxy to value from PropagateToDCA to get correct sign
2866 dcaxy = DCA[0];
2867 // --- dcaz on the other hand is unambiguous
2868 //if(dcaz!=(float)DCA[1]&&debug>6) cout<<"hmm... "<<dcaz<<" "<<DCA[1]<<endl;
2869 //if(debug>6) cout<<"GetXYZ is false, filter bit is "<<fbit<<endl;
2870 }
2871
2872 // if(debug>6)
2873 // {
2874 // cout<<"r[0] is "<<r[0]<<" ";
2875 // cout<<"r[1] is "<<r[1]<<" ";
2876 // cout<<"r[2] is "<<r[2]<<endl;
2877 // cout<<"eventX is "<<eventX<<" ";
2878 // cout<<"eventY is "<<eventY<<" ";
2879 // cout<<"eventZ is "<<eventZ<<endl;
2880 // cout<<"Tdcaxy is "<<Tdcaxy<<" and dcaxy is "<<dcaxy<<" and DCACOV xy is "<<DCA[0]<<endl;
2881 // cout<<"Tdcaz is "<<Tdcaz<<" and dcaz is "<<dcaz<<" and DCACOV z is "<<DCA[1]<<endl;
2882 // }
2883
2884 if((fbit==128||fbit==272)&&!dcaflag)
2885 {
2886 dcaxy = Tdcaxy;
2887 dcaz = Tdcaz;
2888 // if(debug>4)
2889 // {
2890 // cout<<"GetXYZ false but should be true, flipping values"<<endl;
2891 // cout<<"Tdcaxy is "<<Tdcaxy<<" and dcaxy is "<<dcaxy<<" and DCACOV xy is "<<DCA[0]<<endl;
2892 // cout<<"Tdcaz is "<<Tdcaz<<" and dcaz is "<<dcaz<<" and DCACOV z is "<<DCA[1]<<endl;
2893 // }
2894 }
2895
2896 // --- track cuts
2897 if(ncls<nclscut) continue;
2898 if(pt<ptmin||pt>ptmax) continue;
2899 if(fabs(eta)>outeta||fabs(eta)<excleta) continue;
2900 if(dodcacuts&&(fabs(dcaxy)>dcacutxy||fabs(dcaz)>dcacutz)) continue;// problems depending on filter bit
2901
2902 bool pos = charge>0.0;
2903 bool neg = charge<0.0;
2904
2905
2906 if(gRandom) delete gRandom;
2907 gRandom = new TRandom3(0);
2908 gRandom->SetSeed(0);
2909 double rand = gRandom->Rndm();
2910 bool sub1 = false;
2911 bool sub2 = false;
2912 if(rand<0.5) sub1 = true;
2913 else sub2 = true;
2914
2915 if(sub1)
2916 {
2917 ntrkA1++;
2918 if(charge>0) ntrkposA1++;
2919 if(charge<0) ntrknegA1++;
2920 }
2921 if(sub2)
2922 {
2923 ntrkA2++;
2924 if(charge>0) ntrkposA2++;
2925 if(charge<0) ntrknegA2++;
2926 }
2927
2928 if(sub1) // random selection
2929 {
2930 for(int i=1; i<hmax;i++)
2931 {
2932 if(charge>0)
2933 {
2934 a1Xp[i] += cos(i*phi);
2935 a1Yp[i] += sin(i*phi);
2936 }
2937 if(charge<0)
2938 {
2939 a1Xn[i] += cos(i*phi);
2940 a1Yn[i] += sin(i*phi);
2941 }
2942 }
2943 // assign values to array element 0
2944 if(charge>0)
2945 {
2946 a1Xp[0] += 1.0;
2947 a1Yp[0] += pt;
2948 }
2949 if(charge<0)
2950 {
2951 a1Xn[0] += 1.0;
2952 a1Yn[0] += pt;
2953 }
2954 }
2955
2956 if(sub2) // random selection
2957 {
2958 for(int i=1; i<hmax;i++)
2959 {
2960 if(charge>0)
2961 {
2962 a2Xp[i] += cos(i*phi);
2963 a2Yp[i] += sin(i*phi);
2964 }
2965 if(charge<0)
2966 {
2967 a2Xn[i] += cos(i*phi);
2968 a2Yn[i] += sin(i*phi);
2969 }
2970 }
2971 // assign values to array element 0
2972 if(charge>0)
2973 {
2974 a2Xp[0] += 1.0;
2975 a2Yp[0] += pt;
2976 }
2977 if(charge<0)
2978 {
2979 a2Xn[0] += 1.0;
2980 a2Yn[0] += pt;
2981 }
2982 }
2983
2984
2985
2986 // -------------------------------------------------------
2987 // --- now calculate differential cumulants for each track
2988 // -------------------------------------------------------
2989 float trkXpl[hmax], trkYpl[hmax], trkXpr[hmax], trkYpr[hmax];
2990 float trkXnl[hmax], trkYnl[hmax], trkXnr[hmax], trkYnr[hmax];
2991 for(int i=0; i<hmax;i++)
2992 {
2993 trkXpl[i] = 0.0;
2994 trkYpl[i] = 0.0;
2995 trkXpr[i] = 0.0;
2996 trkYpr[i] = 0.0;
2997 trkXnl[i] = 0.0;
2998 trkYnl[i] = 0.0;
2999 trkXnr[i] = 0.0;
3000 trkYnr[i] = 0.0;
3001 }
3002 // very similar to above, except += changed to =
3003 if(pt>ptmin&&pt<ptmax)
3004 {
3005 if(eta>-outeta&&eta<-excleta)
3006 {
3007 // assign values to array elements 1..8
3008 for(int i=1; i<hmax;i++)
3009 {
3010 if(charge>0)
3011 {
3012 trkXpl[i] = cos(i*phi);
3013 trkYpl[i] = sin(i*phi);
3014 }
3015 if(charge<0)
3016 {
3017 trkXnl[i] = cos(i*phi);
3018 trkYnl[i] = sin(i*phi);
3019 }
3020 }
3021 // assign values to array element 0
3022 if(charge>0)
3023 {
3024 trkXpl[0] = 1.0;
3025 trkYpl[0] = pt;
3026 }
3027 if(charge<0)
3028 {
3029 trkXnl[0] = 1.0;
3030 trkYnl[0] = pt;
3031 }
3032 }
3033 //----------------------------
3034 if(eta>excleta&&eta<outeta) // right half of TPC
3035 {
3036 for(int i=1; i<hmax; i++)
3037 {
3038 if(charge>0)
3039 {
3040 trkXpr[i] = cos(i*phi);
3041 trkYpr[i] = sin(i*phi);
3042 }
3043 if(charge<0)
3044 {
3045 trkXnr[i] = cos(i*phi);
3046 trkYnr[i] = sin(i*phi);
3047 }
3048 }
3049 if(charge>0)
3050 {
3051 trkXpr[0] = 1.0;
3052 trkYpr[0] = pt;
3053 }
3054 if(charge<0)
3055 {
3056 trkXnr[0] = 1.0;
3057 trkYnr[0] = pt;
3058 }
3059 } // end right half of TPC
3060 } // end pT selection
3061
3062 // COME BACK HERE
3063 float trkX[9], trkY[9];
3064 float trkXp[9], trkYp[9]; // pos
3065 float trkXn[9], trkYn[9]; // neg
3066 for(int i=0; i<6; i++)
3067 {
3068 trkX[i]=trkXpl[i]+trkXnl[i]+trkXpr[i]+trkXnr[i];
3069 trkY[i]=trkYpl[i]+trkYnl[i]+trkYpr[i]+trkYnr[i];
3070 // pos
3071 trkXp[i]=trkXpl[i]+trkXpr[i];
3072 trkYp[i]=trkYpl[i]+trkYpr[i];
3073 // neg
3074 trkXn[i]=trkXnl[i]+trkXnr[i];
3075 trkYn[i]=trkYnl[i]+trkYnr[i];
3076 }
3077
3078 // sanity checks to prevent division by zero and other problems
3079 if(trkX[0]<1) continue;
3080 if(tpcX[0]<2) continue;
3081 if(tpcXp[0]<2) continue;
3082 if(tpcXn[0]<2) continue;
3083
3084 // --- calculate differential cumulants from Q-vector components
3085 // --- reference flow is from whole TPC
3086 float diffq22ev = ((trkX[2]*tpcX[2])+(trkY[2]*tpcY[2])-1)/(tpcX[0]-1);
3087 float diffq22Pev = ((trkXp[2]*tpcX[2])+(trkYp[2]*tpcY[2])-1)/(tpcX[0]-1);
3088 float diffq22Nev = ((trkXn[2]*tpcX[2])+(trkYn[2]*tpcY[2])-1)/(tpcX[0]-1);
3089 // --- third harmonic
3090 float diffq32ev = ((trkX[3]*tpcX[3])+(trkY[3]*tpcY[3])-1)/(tpcX[0]-1);
3091 float diffq32Pev = ((trkXp[3]*tpcX[3])+(trkYp[3]*tpcY[3])-1)/(tpcX[0]-1);
3092 float diffq32Nev = ((trkXn[3]*tpcX[3])+(trkYn[3]*tpcY[3])-1)/(tpcX[0]-1);
3093 // --- fourth harmonic
3094 float diffq42ev = ((trkX[4]*tpcX[4])+(trkY[4]*tpcY[4])-1)/(tpcX[0]-1);
3095 float diffq42Pev = ((trkXp[4]*tpcX[4])+(trkYp[4]*tpcY[4])-1)/(tpcX[0]-1);
3096 float diffq42Nev = ((trkXn[4]*tpcX[4])+(trkYn[4]*tpcY[4])-1)/(tpcX[0]-1);
3097
3098 // --- differential cumulants test with differing RPs
3099 // float Xdiffq22Pev = ((trkXp[2]*tpcX[2])+(trkYp[2]*tpcY[2])-1)/(tpcX[0]-1);
3100 // float Xdiffq22Nev = ((trkXn[2]*tpcX[2])+(trkYn[2]*tpcY[2])-1)/(tpcX[0]-1);
3101 // float Ydiffq22Pev = ((trkXp[2]*tpcXp[2])+(trkYp[2]*tpcYp[2])-1)/(tpcXp[0]-1);
3102 // float Ydiffq22Nev = ((trkXn[2]*tpcXn[2])+(trkYn[2]*tpcYn[2])-1)/(tpcXn[0]-1);
3103 // float Zdiffq22Pev = ((trkXp[2]*tpcXn[2])+(trkYp[2]*tpcYn[2]))/(tpcXn[0]);
3104 // float Zdiffq22Nev = ((trkXn[2]*tpcXp[2])+(trkYn[2]*tpcYp[2]))/(tpcXp[0]);
3105
3106 // --- histograms to check different reference flow
3107 // if(pos) h_Xdiffq22_ptP->Fill(pt,Xdiffq22Pev);
3108 // if(neg) h_Xdiffq22_ptN->Fill(pt,Xdiffq22Nev);
3109 // if(pos) h_Ydiffq22_ptP->Fill(pt,Ydiffq22Pev);
3110 // if(neg) h_Ydiffq22_ptN->Fill(pt,Ydiffq22Nev);
3111 // if(pos) h_Zdiffq22_ptP->Fill(pt,Zdiffq22Pev);
3112 // if(neg) h_Zdiffq22_ptN->Fill(pt,Zdiffq22Nev);
3113
3114 if(cent>=centlo&&cent<centhi)
3115 {
3116 // --- differential cumulant
3117 h_diffq22_pt->Fill(pt,diffq22ev);
3118 h_diffq22_eta->Fill(eta,diffq22ev);
3119 // --- differential cumulant weighted with qasymm
3120 h_diffAq22_pt->Fill(pt,diffq22ev*qasymm);
3121 h_diffAq22_eta->Fill(eta,diffq22ev*qasymm);
3122 if(pos)
3123 {
3124 h_diffq22_ptP->Fill(pt,diffq22Pev);
3125 h_diffq22_etaP->Fill(eta,diffq22Pev);
3126 h_diffAq22_ptP->Fill(pt,diffq22Pev*qasymm);
3127 h_diffAq22_etaP->Fill(eta,diffq22Pev*qasymm);
3128 }
3129 if(neg)
3130 {
3131 h_diffq22_ptN->Fill(pt,diffq22Nev);
3132 h_diffq22_etaN->Fill(eta,diffq22Nev);
3133 h_diffAq22_ptN->Fill(pt,diffq22Nev*qasymm);
3134 h_diffAq22_etaN->Fill(eta,diffq22Nev*qasymm);
3135 }
3136 // --- same for third harmonic
3137 h_diffq32_pt->Fill(pt,diffq32ev);
3138 h_diffq32_eta->Fill(eta,diffq32ev);
3139 h_diffAq32_pt->Fill(pt,diffq32ev*qasymm);
3140 h_diffAq32_eta->Fill(eta,diffq32ev*qasymm);
3141 if(pos)
3142 {
3143 h_diffq32_ptP->Fill(pt,diffq32Pev);
3144 h_diffq32_etaP->Fill(eta,diffq32Pev);
3145 h_diffAq32_ptP->Fill(pt,diffq32Pev*qasymm);
3146 h_diffAq32_etaP->Fill(eta,diffq32Pev*qasymm);
3147 }
3148 if(neg)
3149 {
3150 h_diffq32_ptN->Fill(pt,diffq32Nev);
3151 h_diffq32_etaN->Fill(eta,diffq32Nev);
3152 h_diffAq32_ptN->Fill(pt,diffq32Nev*qasymm);
3153 h_diffAq32_etaN->Fill(eta,diffq32Nev*qasymm);
3154 }
3155 // --- same for fourth harmonic
3156 h_diffq42_pt->Fill(pt,diffq42ev);
3157 h_diffq42_eta->Fill(eta,diffq42ev);
3158 h_diffAq42_pt->Fill(pt,diffq42ev*qasymm);
3159 h_diffAq42_eta->Fill(eta,diffq42ev*qasymm);
3160 if(pos)
3161 {
3162 h_diffq42_ptP->Fill(pt,diffq42Pev);
3163 h_diffq42_etaP->Fill(eta,diffq42Pev);
3164 h_diffAq42_ptP->Fill(pt,diffq42Pev*qasymm);
3165 h_diffAq42_etaP->Fill(eta,diffq42Pev*qasymm);
3166 }
3167 if(neg)
3168 {
3169 h_diffq42_ptN->Fill(pt,diffq42Nev);
3170 h_diffq42_etaN->Fill(eta,diffq42Nev);
3171 h_diffAq42_ptN->Fill(pt,diffq42Nev*qasymm);
3172 h_diffAq42_etaN->Fill(eta,diffq42Nev*qasymm);
3173 }
3174
3175 // selection to do nested track loop
3176 if(donested)
3177 {
3178 for(int i=0; i<d_ntrk; i++)
3179 {
3180 // ----------------------------------------------------------------------------------------
3181 // --- want delta eta and delta pt dependence of differential cumulant weighted with charge
3182 // --- need nested track loop to get eta1 and eta3
3183 // ----------------------------------------------------------------------------------------
3184 // make sure eta is in range, -99 should be autorejected by this
3185 //if(eta3[i]==-99) continue;
3186 if(eta3[i]<-outeta||eta3[i]>outeta)
3187 {
3188 //cout<<"eta3 out of range!!! "<<eta3[i]<<endl;
3189 continue;
3190 }
3191 if(eta<-outeta||eta>outeta)
3192 {
3193 cout<<"eta out of range!!! "<<eta<<endl;
3194 continue;
3195 }
3196 // charge==0 should be rejected by the above cut
3197 if(charge3[i]==0) {cout<<"WTF!"<<endl;continue;}
3198 // VERY IMPORTANT remove auto-correlations with 3rd particle
3199 if(i==itrk) continue;
3200
3201 float DETA = eta-eta3[i];
3202 float DPHI = phi-phi3[i];
3203 float DPT = pt-pt3[i];
3204
3205 // magnetic field sign
3206 //bSign = (eventMain->GetMagneticField() > 0) ? 1 : -1;
3207 float bSign = 1;
3208 if(mag<0) bSign = -1;
3209
3210 // double deta = firstEta - secondEta[j];
3211 // double dphi = firstPhi - secondPhi[j];
3212 float deta = DETA;
3213 float dphi = DPHI;
3214 float fHBTCutValue = 0.02;
3215 // optimization
3216 if(dopaircut&&(fabs(deta)<fHBTCutValue*2.5*3)) //fHBTCutValue = 0.02 [default for dphicorrelations]
3217 {
3218 if(debug>4)
3219 {
3220 cout<<"inside deta cut, evaluating"<<endl;
3221 cout<<"deta is "<<deta<<" and dphi is "<<dphi<<endl;
3222 }
3223 // phi in rad
3224 // float phi1rad = firstPhi;
3225 // float phi2rad = secondPhi[j];
3226
3227 // check first boundaries to see if is worth to loop and find the minimum
3228 // float dphistar1 = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, 0.8, bSign);
3229 // float dphistar2 = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, 2.5, bSign);
3230 float dphistar1 = GetDPhiStar(phi, pt, charge, phi3[i], pt3[i], charge3[i], 0.8, bSign);
3231 float dphistar2 = GetDPhiStar(phi, pt, charge, phi3[i], pt3[i], charge3[i], 2.5, bSign);
3232
3233 const float kLimit = fHBTCutValue * 3;
3234
3235 float dphistarminabs = 1e5;
3236 float dphistarmin = 1e5;
3237
3238 if(fabs(dphistar1) < kLimit || fabs(dphistar2) < kLimit || dphistar1 * dphistar2 < 0 )
3239 {
3240 if(debug>4)
3241 {
3242 cout<<"inside detadphi cut, evaluating"<<endl;
3243 cout<<"deta is "<<deta<<" and dphi is "<<dphi<<endl;
3244 }
3245 for(double rad=0.8; rad<2.51; rad+=0.01)
3246 {
3247 //float dphistar = GetDPhiStar(phi1rad, firstPt, charge1, phi2rad, secondPt[j], charge2, rad, bSign);
3248 float dphistar = GetDPhiStar(phi, pt, charge, phi3[i], pt3[i], charge3[i], 0.8, bSign);
3249 float dphistarabs = fabs(dphistar);
3250 if(dphistarabs < dphistarminabs)
3251 {
3252 dphistarmin = dphistar;
3253 dphistarminabs = dphistarabs;
3254 }
3255 }
3256 if(dphistarminabs < fHBTCutValue && fabs(deta) < fHBTCutValue)
3257 {
3258 if(debug>4)
3259 printf("HBT: Removed track pair %d %d with [[%f %f]] %f %f %f | %f %f %d %f %f %d %f",
3260 itrk, i, deta, dphi, dphistarminabs, dphistar1, dphistar2, phi, pt, charge, phi3[i], pt3[i], charge3[i], bSign);
3261 continue;
3262 }
3263 } // loop over angles
3264 } // selection on pairs close together
3265
3266
3267
3268
3269
3270 h_AT_X_deta->Fill(DETA,charge3[i]);
3271 if(pos) h_AT_X_detaP->Fill(DETA,charge3[i]);
3272 if(neg) h_AT_X_detaN->Fill(DETA,charge3[i]);
3273 h_AT_X_dpt->Fill(DPT,charge3[i]);
3274 if(pos) h_AT_X_dptP->Fill(DPT,charge3[i]);
3275 if(neg) h_AT_X_dptN->Fill(DPT,charge3[i]);
3276 // --- particles 1, 2, and 3 in whole TPC
3277 if(pos) h_diffq22_X_detaP->Fill(DETA,diffq22Pev);
3278 if(neg) h_diffq22_X_detaN->Fill(DETA,diffq22Nev);
3279 if(pos) h_diffq32_X_detaP->Fill(DETA,diffq32Pev);
3280 if(neg) h_diffq32_X_detaN->Fill(DETA,diffq32Nev);
3281 if(pos) h_diffq42_X_detaP->Fill(DETA,diffq42Pev);
3282 if(neg) h_diffq42_X_detaN->Fill(DETA,diffq42Nev);
3283 if(pos) h_ATdiffq22_X_detaP->Fill(DETA,diffq22Pev*charge3[i]);
3284 if(neg) h_ATdiffq22_X_detaN->Fill(DETA,diffq22Nev*charge3[i]);
3285 if(pos) h_ATdiffq32_X_detaP->Fill(DETA,diffq32Pev*charge3[i]);
3286 if(neg) h_ATdiffq32_X_detaN->Fill(DETA,diffq32Nev*charge3[i]);
3287 if(pos) h_ATdiffq42_X_detaP->Fill(DETA,diffq42Pev*charge3[i]);
3288 if(neg) h_ATdiffq42_X_detaN->Fill(DETA,diffq42Nev*charge3[i]);
3289 if(pos) h_diffq22_X_dptP->Fill(DPT,diffq22Pev);
3290 if(neg) h_diffq22_X_dptN->Fill(DPT,diffq22Nev);
3291 if(pos) h_diffq32_X_dptP->Fill(DPT,diffq32Pev);
3292 if(neg) h_diffq32_X_dptN->Fill(DPT,diffq32Nev);
3293 if(pos) h_diffq42_X_dptP->Fill(DPT,diffq42Pev);
3294 if(neg) h_diffq42_X_dptN->Fill(DPT,diffq42Nev);
3295 if(pos) h_ATdiffq22_X_dptP->Fill(DPT,diffq22Pev*charge3[i]);
3296 if(neg) h_ATdiffq22_X_dptN->Fill(DPT,diffq22Nev*charge3[i]);
3297 if(pos) h_ATdiffq32_X_dptP->Fill(DPT,diffq32Pev*charge3[i]);
3298 if(neg) h_ATdiffq32_X_dptN->Fill(DPT,diffq32Nev*charge3[i]);
3299 if(pos) h_ATdiffq42_X_dptP->Fill(DPT,diffq42Pev*charge3[i]);
3300 if(neg) h_ATdiffq42_X_dptN->Fill(DPT,diffq42Nev*charge3[i]);
3301
3302 float meanchargeeta = 0;
3303 int index = int((eta3[i]+1.0)*50);
3304 if(index<0||index>100)
3305 {
3306 cout<<"out of index!!!"<<endl;
3307 continue; // safety, should be redundant
3308 }
3309 //cout<<"eta is "<<eta3[i]<<" index is "<<index<<endl;
3310 if(mag<-4.0)
3311 {
3312 if(fbit==1) meanchargeeta = MCEF1_fb1[index];
3313 if(fbit==128) meanchargeeta = MCEF1_fb128[index];
3314 if(fbit==272) meanchargeeta = MCEF1_fb272[index];
3315 //if(fbit==1&&dodcacuts) meanchargeeta = MCEF1_dca[index];
3316 }
3317 if(mag>4.0)
3318 {
3319 if(fbit==1) meanchargeeta = MCEF3_fb1[index];
3320 if(fbit==128) meanchargeeta = MCEF3_fb128[index];
3321 if(fbit==272) meanchargeeta = MCEF3_fb272[index];
3322 //if(fbit==1&&dodcacuts) meanchargeeta = MCEF3_dca[index];
3323 }
3324 float charge3subeta = charge3[i] - meanchargeeta;
3325
3326 h_AT_S_deta->Fill(DETA,charge3subeta);
3327 if(pos) h_AT_S_detaP->Fill(DETA,charge3subeta);
3328 if(neg) h_AT_S_detaN->Fill(DETA,charge3subeta);
3329 // --- particles 1, 2, and 3 in whole TPC
3330 if(pos) h_diffq22_S_detaP->Fill(DETA,diffq22Pev);
3331 if(neg) h_diffq22_S_detaN->Fill(DETA,diffq22Nev);
3332 if(pos) h_diffq32_S_detaP->Fill(DETA,diffq32Pev);
3333 if(neg) h_diffq32_S_detaN->Fill(DETA,diffq32Nev);
3334 if(pos) h_diffq42_S_detaP->Fill(DETA,diffq42Pev);
3335 if(neg) h_diffq42_S_detaN->Fill(DETA,diffq42Nev);
3336 if(pos) h_ATdiffq22_S_detaP->Fill(DETA,diffq22Pev*charge3subeta);
3337 if(neg) h_ATdiffq22_S_detaN->Fill(DETA,diffq22Nev*charge3subeta);
3338 if(pos) h_ATdiffq32_S_detaP->Fill(DETA,diffq32Pev*charge3subeta);
3339 if(neg) h_ATdiffq32_S_detaN->Fill(DETA,diffq32Nev*charge3subeta);
3340 if(pos) h_ATdiffq42_S_detaP->Fill(DETA,diffq42Pev*charge3subeta);
3341 if(neg) h_ATdiffq42_S_detaN->Fill(DETA,diffq42Nev*charge3subeta);
3342
3343 // --- now field selection for delta eta stuff
3344 if(mag<-4.0)
3345 {
3346 h_AT_X_deta_F1->Fill(DETA,charge3[i]);
3347 if(pos) h_AT_X_detaP_F1->Fill(DETA,charge3[i]);
3348 if(neg) h_AT_X_detaN_F1->Fill(DETA,charge3[i]);
3349 if(pos) h_diffq22_X_detaP_F1->Fill(DETA,diffq22Pev);
3350 if(neg) h_diffq22_X_detaN_F1->Fill(DETA,diffq22Nev);
3351 if(pos) h_diffq32_X_detaP_F1->Fill(DETA,diffq32Pev);
3352 if(neg) h_diffq32_X_detaN_F1->Fill(DETA,diffq32Nev);
3353 if(pos) h_diffq42_X_detaP_F1->Fill(DETA,diffq42Pev);
3354 if(neg) h_diffq42_X_detaN_F1->Fill(DETA,diffq42Nev);
3355 if(pos) h_ATdiffq22_X_detaP_F1->Fill(DETA,diffq22Pev*charge3[i]);
3356 if(neg) h_ATdiffq22_X_detaN_F1->Fill(DETA,diffq22Nev*charge3[i]);
3357 if(pos) h_ATdiffq32_X_detaP_F1->Fill(DETA,diffq32Pev*charge3[i]);
3358 if(neg) h_ATdiffq32_X_detaN_F1->Fill(DETA,diffq32Nev*charge3[i]);
3359 if(pos) h_ATdiffq42_X_detaP_F1->Fill(DETA,diffq42Pev*charge3[i]);
3360 if(neg) h_ATdiffq42_X_detaN_F1->Fill(DETA,diffq42Nev*charge3[i]);
3361 }
3362 if(mag>4.0)
3363 {
3364 h_AT_X_deta_F3->Fill(DETA,charge3[i]);
3365 if(pos) h_AT_X_detaP_F3->Fill(DETA,charge3[i]);
3366 if(neg) h_AT_X_detaN_F3->Fill(DETA,charge3[i]);
3367 if(pos) h_diffq22_X_detaP_F3->Fill(DETA,diffq22Pev);
3368 if(neg) h_diffq22_X_detaN_F3->Fill(DETA,diffq22Nev);
3369 if(pos) h_diffq32_X_detaP_F3->Fill(DETA,diffq32Pev);
3370 if(neg) h_diffq32_X_detaN_F3->Fill(DETA,diffq32Nev);
3371 if(pos) h_diffq42_X_detaP_F3->Fill(DETA,diffq42Pev);
3372 if(neg) h_diffq42_X_detaN_F3->Fill(DETA,diffq42Nev);
3373 if(pos) h_ATdiffq22_X_detaP_F3->Fill(DETA,diffq22Pev*charge3[i]);
3374 if(neg) h_ATdiffq22_X_detaN_F3->Fill(DETA,diffq22Nev*charge3[i]);
3375 if(pos) h_ATdiffq32_X_detaP_F3->Fill(DETA,diffq32Pev*charge3[i]);
3376 if(neg) h_ATdiffq32_X_detaN_F3->Fill(DETA,diffq32Nev*charge3[i]);
3377 if(pos) h_ATdiffq42_X_detaP_F3->Fill(DETA,diffq42Pev*charge3[i]);
3378 if(neg) h_ATdiffq42_X_detaN_F3->Fill(DETA,diffq42Nev*charge3[i]);
3379 }
3380
3381 } // end nested track loop
3382
3383 } // end if statement to check for nested track loop
3384
3385 } // centrality selection for differential cumulants
3386
3387 } // end of second track loop
3388
3389 float qasymA1 = (float)(ntrkposA1-ntrknegA1)/(ntrkA1);
3390 float qasymA2 = (float)(ntrkposA2-ntrknegA2)/(ntrkA2);
3391
3392 h_rA_cent->Fill(cent,qasymA1);
3393 h_r1A_cent->Fill(cent,qasymA1);
3394 h_r2A_cent->Fill(cent,qasymA2);
3395
3396 h2_rA_cent->Fill(cent,qasymA1);
3397 h2_r1A_cent->Fill(cent,qasymA1);
3398 h2_r2A_cent->Fill(cent,qasymA2);
3399
3400 // random subevents 1
3401 float a1XX[9], a1YY[9], a1QQ[9];//, a1qq[9];
3402 float a1XXp[9], a1YYp[9], a1QQp[9];//, a1qqp[9]; // pos
3403 float a1XXn[9], a1YYn[9], a1QQn[9];//, a1qqn[9]; // neg
3404 for(int i=0; i<6; i++)
3405 {
3406 a1XX[i]=a1Xp[i]+a1Xn[i];
3407 a1YY[i]=a1Yp[i]+a1Yn[i];
3408 a1QQ[i]=a1XX[i]*a1XX[i]+a1YY[i]*a1YY[i];
3409 //a1qq[i]=sqrt(a1QQ[i]/a1XX[0]);
3410 // pos
3411 a1XXp[i]=a1Xp[i];
3412 a1YYp[i]=a1Yp[i];
3413 a1QQp[i]=a1XXp[i]*a1XXp[i]+a1YYp[i]*a1YYp[i];
3414 //a1qqp[i]=sqrt(a1QQp[i]/a1XXp[0]);
3415 // neg
3416 a1XXn[i]=a1Xn[i];
3417 a1YYn[i]=a1Yn[i];
3418 a1QQn[i]=a1XXn[i]*a1XXn[i]+a1YYn[i]*a1YYn[i];
3419 //a1qqn[i]=sqrt(a1QQn[i]/a1XXn[0]);
3420 }
3421 // random subevents 2
3422 float a2XX[9], a2YY[9], a2QQ[9];//, a2qq[9];
3423 float a2XXp[9], a2YYp[9], a2QQp[9];//, a2qqp[9]; // pos
3424 float a2XXn[9], a2YYn[9], a2QQn[9];//, a2qqn[9]; // neg
3425 for(int i=0; i<6; i++)
3426 {
3427 a2XX[i]=a2Xp[i]+a2Xn[i];
3428 a2YY[i]=a2Yp[i]+a2Yn[i];
3429 a2QQ[i]=a2XX[i]*a2XX[i]+a2YY[i]*a2YY[i];
3430 //a2qq[i]=sqrt(a2QQ[i]/a2XX[0]);
3431 // pos
3432 a2XXp[i]=a2Xp[i];
3433 a2YYp[i]=a2Yp[i];
3434 a2QQp[i]=a2XXp[i]*a2XXp[i]+a2YYp[i]*a2YYp[i];
3435 //a2qqp[i]=sqrt(a2QQp[i]/a2XXp[0]);
3436 // neg
3437 a2XXn[i]=a2Xn[i];
3438 a2YYn[i]=a2Yn[i];
3439 a2QQn[i]=a2XXn[i]*a2XXn[i]+a2YYn[i]*a2YYn[i];
3440 //a2qqn[i]=sqrt(a1QQn[i]/a1XXn[0]);
3441 }
3442
3443 float a1M = a1XX[0];
3444 float a1W_2 = a1M*(a1M-(a1M/M));
3445 float a1Mp = a1XXp[0];
3446 float a1Wp_2 = a1Mp*(a1Mp-(a1Mp/Mp));
3447 float a1Mn = a1XXn[0];
3448 float a1Wn_2 = a1Mn*(a1Mn-(a1Mn/Mn));
3449
3450 float a2M = a2XX[0];
3451 float a2W_2 = a2M*(a2M-(a2M/M));
3452 float a2Mp = a2XXp[0];
3453 float a2Wp_2 = a2Mp*(a2Mp-(a2Mp/Mp));
3454 float a2Mn = a2XXn[0];
3455 float a2Wn_2 = a2Mn*(a2Mn-(a2Mn/Mn));
3456
3457 float a1q22ev = (a1QQ[2]-a1M)/a1W_2;
3458 float a1q22Pev = (a1QQp[2]-a1Mp)/a1Wp_2; // pos
3459 float a1q22Nev = (a1QQn[2]-a1Mn)/a1Wn_2; // neg
3460
3461 float a2q22ev = (a2QQ[2]-a2M)/a2W_2;
3462 float a2q22Pev = (a2QQp[2]-a2Mp)/a2Wp_2; // pos
3463 float a2q22Nev = (a2QQn[2]-a2Mn)/a2Wn_2; // neg
3464
3465 h_a1q22_cent->Fill(cent,a1q22ev);
3466 h_a1q22_centP->Fill(cent,a1q22Pev);
3467 h_a1q22_centN->Fill(cent,a1q22Nev);
3468
3469 h_a2q22_cent->Fill(cent,a2q22ev);
3470 h_a2q22_centP->Fill(cent,a2q22Pev);
3471 h_a2q22_centN->Fill(cent,a2q22Nev);
3472
3473 h_rAq22_X1_cent->Fill(cent,a1q22ev*qasymA1);
3474 h_rAq22_X1_centP->Fill(cent,a1q22Pev*qasymA1);
3475 h_rAq22_X1_centN->Fill(cent,a1q22Nev*qasymA1);
3476 h_rAq22_X2_cent->Fill(cent,a2q22ev*qasymA2);
3477 h_rAq22_X2_centP->Fill(cent,a2q22Pev*qasymA2);
3478 h_rAq22_X2_centN->Fill(cent,a2q22Nev*qasymA2);
3479 h_rAq22_X3_cent->Fill(cent,a1q22ev*qasymA2);
3480 h_rAq22_X3_centP->Fill(cent,a1q22Pev*qasymA2);
3481 h_rAq22_X3_centN->Fill(cent,a1q22Nev*qasymA2);
3482 h_rAq22_X4_cent->Fill(cent,a2q22ev*qasymA1);
3483 h_rAq22_X4_centP->Fill(cent,a2q22Pev*qasymA1);
3484 h_rAq22_X4_centN->Fill(cent,a2q22Nev*qasymA1);
3485
3486
3487
3488 // ------------------------------------ //
3489 // --- send data to the output list --- //
3490 // ------------------------------------ //
3491
3492 PostData(1,fOutputList);
3493
3494 // ------------------------- //
3495 // --- end of event loop --- //
3496 // ------------------------- //
3497
3498}
3499
3500void AliAnalysisTaskCMEv2A::SetParameters()
3501{
3502
3503 cout<<"Now setting parameters!"<<endl;
3504
3505 // --- these are the TPC efficiencies
3506 effTPC[0] = 0;//pt=0.05
3507 effTPC[1] = 0;//pt=0.15
3508 effTPC[2] = 0.895961;//pt=0.25
3509 effTPC[3] = 0.92504;//pt=0.35
3510 effTPC[4] = 0.944326;//pt=0.45
3511 effTPC[5] = 0.948261;//pt=0.55
3512 effTPC[6] = 0.948522;//pt=0.65
3513 effTPC[7] = 0.946941;//pt=0.75
3514 effTPC[8] = 0.948537;//pt=0.85
3515 effTPC[9] = 0.95194;//pt=0.95
3516 effTPC[10] = 0.95228;//pt=1.05
3517 effTPC[11] = 0.954775;//pt=1.15
3518 effTPC[12] = 0.956055;//pt=1.25
3519 effTPC[13] = 0.955758;//pt=1.35
3520 effTPC[14] = 0.951268;//pt=1.45
3521 effTPC[15] = 0.952081;//pt=1.55
3522 effTPC[16] = 0.947991;//pt=1.65
3523 effTPC[17] = 0.943755;//pt=1.75
3524 effTPC[18] = 0.938106;//pt=1.85
3525 effTPC[19] = 0.932583;//pt=1.95
3526 effTPC[20] = 0.927206;//pt=2.05
3527 effTPC[21] = 0.921111;//pt=2.15
3528 effTPC[22] = 0.913127;//pt=2.25
3529 effTPC[23] = 0.910508;//pt=2.35
3530 effTPC[24] = 0.905838;//pt=2.45
3531 effTPC[25] = 0.90243;//pt=2.55
3532 effTPC[26] = 0.897893;//pt=2.65
3533 effTPC[27] = 0.895459;//pt=2.75
3534 effTPC[28] = 0.89308;//pt=2.85
3535 effTPC[29] = 0.895013;//pt=2.95
3536 effTPC[30] = 0.890872;//pt=3.05
3537 effTPC[31] = 0.892046;//pt=3.15
3538 effTPC[32] = 0.890133;//pt=3.25
3539 effTPC[33] = 0.894786;//pt=3.35
3540 effTPC[34] = 0.891886;//pt=3.45
3541 effTPC[35] = 0.893649;//pt=3.55
3542 effTPC[36] = 0.891418;//pt=3.65
3543 effTPC[37] = 0.897408;//pt=3.75
3544 effTPC[38] = 0.890273;//pt=3.85
3545 effTPC[39] = 0.896526;//pt=3.95
3546 effTPC[40] = 0.892236;//pt=4.05
3547 effTPC[41] = 0.895075;//pt=4.15
3548 effTPC[42] = 0.895159;//pt=4.25
3549 effTPC[43] = 0.889562;//pt=4.35
3550 effTPC[44] = 0.890777;//pt=4.45
3551 effTPC[45] = 0.892925;//pt=4.55
3552 effTPC[46] = 0.890358;//pt=4.65
3553 effTPC[47] = 0.897014;//pt=4.75
3554 effTPC[48] = 0.886445;//pt=4.85
3555 effTPC[49] = 0.890928;//pt=4.95
3556
3557 float max = effTPC[0];
3558 int maxbin = 0;
3559 for(int i=0; i<49; i++)
3560 {
3561 if(effTPC[i] > max)
3562 {
3563 max = effTPC[i];
3564 maxbin = i;
3565 }
3566 }
3567 cout<<"max for TPC is "<<max<<endl;
3568 cout<<"maxbin for TPC is "<<maxbin<<endl;
3569 for(int i=0; i<50; i++)
3570 {
3571 effTPC[i] /= max;
3572 }
3573
3574 effHYB[0] = 0;//pt=0.05
3575 effHYB[1] = 0;//pt=0.15
3576 effHYB[2] = 0.797902;//pt=0.25
3577 effHYB[3] = 0.842964;//pt=0.35
3578 effHYB[4] = 0.872516;//pt=0.45
3579 effHYB[5] = 0.883225;//pt=0.55
3580 effHYB[6] = 0.886472;//pt=0.65
3581 effHYB[7] = 0.886011;//pt=0.75
3582 effHYB[8] = 0.883906;//pt=0.85
3583 effHYB[9] = 0.882904;//pt=0.95
3584 effHYB[10] = 0.881162;//pt=1.05
3585 effHYB[11] = 0.883611;//pt=1.15
3586 effHYB[12] = 0.884113;//pt=1.25
3587 effHYB[13] = 0.883676;//pt=1.35
3588 effHYB[14] = 0.875165;//pt=1.45
3589 effHYB[15] = 0.87379;//pt=1.55
3590 effHYB[16] = 0.871452;//pt=1.65
3591 effHYB[17] = 0.868833;//pt=1.75
3592 effHYB[18] = 0.863982;//pt=1.85
3593 effHYB[19] = 0.858165;//pt=1.95
3594 effHYB[20] = 0.850305;//pt=2.05
3595 effHYB[21] = 0.845429;//pt=2.15
3596 effHYB[22] = 0.839418;//pt=2.25
3597 effHYB[23] = 0.836466;//pt=2.35
3598 effHYB[24] = 0.833551;//pt=2.45
3599 effHYB[25] = 0.831225;//pt=2.55
3600 effHYB[26] = 0.826823;//pt=2.65
3601 effHYB[27] = 0.822582;//pt=2.75
3602 effHYB[28] = 0.82375;//pt=2.85
3603 effHYB[29] = 0.824866;//pt=2.95
3604 effHYB[30] = 0.822234;//pt=3.05
3605 effHYB[31] = 0.825215;//pt=3.15
3606 effHYB[32] = 0.82225;//pt=3.25
3607 effHYB[33] = 0.825158;//pt=3.35
3608 effHYB[34] = 0.822447;//pt=3.45
3609 effHYB[35] = 0.825239;//pt=3.55
3610 effHYB[36] = 0.825619;//pt=3.65
3611 effHYB[37] = 0.831917;//pt=3.75
3612 effHYB[38] = 0.82425;//pt=3.85
3613 effHYB[39] = 0.830073;//pt=3.95
3614 effHYB[40] = 0.82808;//pt=4.05
3615 effHYB[41] = 0.829005;//pt=4.15
3616 effHYB[42] = 0.831138;//pt=4.25
3617 effHYB[43] = 0.827304;//pt=4.35
3618 effHYB[44] = 0.830745;//pt=4.45
3619 effHYB[45] = 0.832802;//pt=4.55
3620 effHYB[46] = 0.828073;//pt=4.65
3621 effHYB[47] = 0.832246;//pt=4.75
3622 effHYB[48] = 0.82791;//pt=4.85
3623 effHYB[49] = 0.831579;//pt=4.95
3624
3625 max = effHYB[0];
3626 maxbin = 0;
3627 for(int i=0; i<49; i++)
3628 {
3629 if(effHYB[i] > max)
3630 {
3631 max = effHYB[i];
3632 maxbin = i;
3633 }
3634 }
3635 cout<<"max for HYB is "<<max<<endl;
3636 cout<<"maxbin for HYB is "<<maxbin<<endl;
3637 for(int i=0; i<50; i++)
3638 {
3639 effHYB[i] /= max;
3640 }
3641
3642 // ------------------------------------
3643 // --- now for mean charge subtractions
3644 // ------------------------------------
3645
3646 // --- filter bit 1 (global tracks)
3647 MCEF1_fb1[0] = 0.00591299; MCEF3_fb1[0] = 0.00383748;//eta is -0.99
3648 MCEF1_fb1[1] = 0.00588992; MCEF3_fb1[1] = 0.003667;//eta is -0.97
3649 MCEF1_fb1[2] = 0.00570083; MCEF3_fb1[2] = 0.00361253;//eta is -0.95
3650 MCEF1_fb1[3] = 0.00560151; MCEF3_fb1[3] = 0.00333053;//eta is -0.93
3651 MCEF1_fb1[4] = 0.00576251; MCEF3_fb1[4] = 0.00340385;//eta is -0.91
3652 MCEF1_fb1[5] = 0.00554351; MCEF3_fb1[5] = 0.00292592;//eta is -0.89
3653 MCEF1_fb1[6] = 0.0056127; MCEF3_fb1[6] = 0.00332485;//eta is -0.87
3654 MCEF1_fb1[7] = 0.00541225; MCEF3_fb1[7] = 0.00331356;//eta is -0.85
3655 MCEF1_fb1[8] = 0.00525774; MCEF3_fb1[8] = 0.00370341;//eta is -0.83
3656 MCEF1_fb1[9] = 0.00535911; MCEF3_fb1[9] = 0.00380439;//eta is -0.81
3657 MCEF1_fb1[10] = 0.00532029; MCEF3_fb1[10] = 0.00418114;//eta is -0.79
3658 MCEF1_fb1[11] = 0.00542025; MCEF3_fb1[11] = 0.00387263;//eta is -0.77
3659 MCEF1_fb1[12] = 0.00520206; MCEF3_fb1[12] = 0.00414251;//eta is -0.75
3660 MCEF1_fb1[13] = 0.00513343; MCEF3_fb1[13] = 0.00400127;//eta is -0.73
3661 MCEF1_fb1[14] = 0.00543818; MCEF3_fb1[14] = 0.00416512;//eta is -0.71
3662 MCEF1_fb1[15] = 0.00553948; MCEF3_fb1[15] = 0.00406533;//eta is -0.69
3663 MCEF1_fb1[16] = 0.00556544; MCEF3_fb1[16] = 0.00428908;//eta is -0.67
3664 MCEF1_fb1[17] = 0.00558886; MCEF3_fb1[17] = 0.00435018;//eta is -0.65
3665 MCEF1_fb1[18] = 0.0057369; MCEF3_fb1[18] = 0.00445061;//eta is -0.63
3666 MCEF1_fb1[19] = 0.00616797; MCEF3_fb1[19] = 0.00443836;//eta is -0.61
3667 MCEF1_fb1[20] = 0.00604426; MCEF3_fb1[20] = 0.00466612;//eta is -0.59
3668 MCEF1_fb1[21] = 0.00581836; MCEF3_fb1[21] = 0.00473329;//eta is -0.57
3669 MCEF1_fb1[22] = 0.00573374; MCEF3_fb1[22] = 0.00429524;//eta is -0.55
3670 MCEF1_fb1[23] = 0.00563837; MCEF3_fb1[23] = 0.0047092;//eta is -0.53
3671 MCEF1_fb1[24] = 0.00628575; MCEF3_fb1[24] = 0.0050631;//eta is -0.51
3672 MCEF1_fb1[25] = 0.00606812; MCEF3_fb1[25] = 0.00480595;//eta is -0.49
3673 MCEF1_fb1[26] = 0.00614111; MCEF3_fb1[26] = 0.00474149;//eta is -0.47
3674 MCEF1_fb1[27] = 0.0060786; MCEF3_fb1[27] = 0.00478341;//eta is -0.45
3675 MCEF1_fb1[28] = 0.00621897; MCEF3_fb1[28] = 0.00499422;//eta is -0.43
3676 MCEF1_fb1[29] = 0.00633416; MCEF3_fb1[29] = 0.00498704;//eta is -0.41
3677 MCEF1_fb1[30] = 0.00638646; MCEF3_fb1[30] = 0.00503349;//eta is -0.39
3678 MCEF1_fb1[31] = 0.00658823; MCEF3_fb1[31] = 0.00485434;//eta is -0.37
3679 MCEF1_fb1[32] = 0.0064771; MCEF3_fb1[32] = 0.00482841;//eta is -0.35
3680 MCEF1_fb1[33] = 0.00694801; MCEF3_fb1[33] = 0.00471101;//eta is -0.33
3681 MCEF1_fb1[34] = 0.00636795; MCEF3_fb1[34] = 0.00491829;//eta is -0.31
3682 MCEF1_fb1[35] = 0.00650891; MCEF3_fb1[35] = 0.00498886;//eta is -0.29
3683 MCEF1_fb1[36] = 0.00653225; MCEF3_fb1[36] = 0.00481843;//eta is -0.27
3684 MCEF1_fb1[37] = 0.00657161; MCEF3_fb1[37] = 0.00511179;//eta is -0.25
3685 MCEF1_fb1[38] = 0.00640258; MCEF3_fb1[38] = 0.00511465;//eta is -0.23
3686 MCEF1_fb1[39] = 0.00626697; MCEF3_fb1[39] = 0.00563881;//eta is -0.21
3687 MCEF1_fb1[40] = 0.00691174; MCEF3_fb1[40] = 0.00534971;//eta is -0.19
3688 MCEF1_fb1[41] = 0.00620562; MCEF3_fb1[41] = 0.00557349;//eta is -0.17
3689 MCEF1_fb1[42] = 0.00649612; MCEF3_fb1[42] = 0.00498584;//eta is -0.15
3690 MCEF1_fb1[43] = 0.0065798; MCEF3_fb1[43] = 0.00542058;//eta is -0.13
3691 MCEF1_fb1[44] = 0.00642303; MCEF3_fb1[44] = 0.00613666;//eta is -0.11
3692 MCEF1_fb1[45] = 0.00653218; MCEF3_fb1[45] = 0.00546797;//eta is -0.09
3693 MCEF1_fb1[46] = 0.0064394; MCEF3_fb1[46] = 0.00555092;//eta is -0.07
3694 MCEF1_fb1[47] = 0.00622478; MCEF3_fb1[47] = 0.00629308;//eta is -0.05
3695 MCEF1_fb1[48] = 0.00639885; MCEF3_fb1[48] = 0.0063706;//eta is -0.03
3696 MCEF1_fb1[49] = 0.00666821; MCEF3_fb1[49] = 0.00631535;//eta is -0.01
3697 MCEF1_fb1[50] = 0.00653648; MCEF3_fb1[50] = 0.00710942;//eta is 0.01
3698 MCEF1_fb1[51] = 0.00638104; MCEF3_fb1[51] = 0.00709238;//eta is 0.03
3699 MCEF1_fb1[52] = 0.0066302; MCEF3_fb1[52] = 0.00752929;//eta is 0.05
3700 MCEF1_fb1[53] = 0.00660631; MCEF3_fb1[53] = 0.00756394;//eta is 0.07
3701 MCEF1_fb1[54] = 0.0068361; MCEF3_fb1[54] = 0.00763778;//eta is 0.09
3702 MCEF1_fb1[55] = 0.00680126; MCEF3_fb1[55] = 0.00766098;//eta is 0.11
3703 MCEF1_fb1[56] = 0.00673586; MCEF3_fb1[56] = 0.0077562;//eta is 0.13
3704 MCEF1_fb1[57] = 0.0066831; MCEF3_fb1[57] = 0.00749726;//eta is 0.15
3705 MCEF1_fb1[58] = 0.00632194; MCEF3_fb1[58] = 0.00814879;//eta is 0.17
3706 MCEF1_fb1[59] = 0.00615406; MCEF3_fb1[59] = 0.00773741;//eta is 0.19
3707 MCEF1_fb1[60] = 0.00640291; MCEF3_fb1[60] = 0.00709483;//eta is 0.21
3708 MCEF1_fb1[61] = 0.00626011; MCEF3_fb1[61] = 0.00748346;//eta is 0.23
3709 MCEF1_fb1[62] = 0.00588407; MCEF3_fb1[62] = 0.00746276;//eta is 0.25
3710 MCEF1_fb1[63] = 0.00615514; MCEF3_fb1[63] = 0.0073329;//eta is 0.27
3711 MCEF1_fb1[64] = 0.00581758; MCEF3_fb1[64] = 0.00749753;//eta is 0.29
3712 MCEF1_fb1[65] = 0.00552913; MCEF3_fb1[65] = 0.00758203;//eta is 0.31
3713 MCEF1_fb1[66] = 0.00618476; MCEF3_fb1[66] = 0.00736015;//eta is 0.33
3714 MCEF1_fb1[67] = 0.00571961; MCEF3_fb1[67] = 0.00717504;//eta is 0.35
3715 MCEF1_fb1[68] = 0.00596315; MCEF3_fb1[68] = 0.00731385;//eta is 0.37
3716 MCEF1_fb1[69] = 0.00617447; MCEF3_fb1[69] = 0.00715388;//eta is 0.39
3717 MCEF1_fb1[70] = 0.00548738; MCEF3_fb1[70] = 0.00691641;//eta is 0.41
3718 MCEF1_fb1[71] = 0.00526914; MCEF3_fb1[71] = 0.00712093;//eta is 0.43
3719 MCEF1_fb1[72] = 0.00549674; MCEF3_fb1[72] = 0.00689862;//eta is 0.45
3720 MCEF1_fb1[73] = 0.00556518; MCEF3_fb1[73] = 0.00701435;//eta is 0.47
3721 MCEF1_fb1[74] = 0.005665; MCEF3_fb1[74] = 0.00676622;//eta is 0.49
3722 MCEF1_fb1[75] = 0.00552991; MCEF3_fb1[75] = 0.00679297;//eta is 0.51
3723 MCEF1_fb1[76] = 0.00569837; MCEF3_fb1[76] = 0.00662394;//eta is 0.53
3724 MCEF1_fb1[77] = 0.00562372; MCEF3_fb1[77] = 0.00638628;//eta is 0.55
3725 MCEF1_fb1[78] = 0.00543001; MCEF3_fb1[78] = 0.00657229;//eta is 0.57
3726 MCEF1_fb1[79] = 0.00524306; MCEF3_fb1[79] = 0.00635786;//eta is 0.59
3727 MCEF1_fb1[80] = 0.00535293; MCEF3_fb1[80] = 0.00641706;//eta is 0.61
3728 MCEF1_fb1[81] = 0.00514873; MCEF3_fb1[81] = 0.00636306;//eta is 0.63
3729 MCEF1_fb1[82] = 0.00510695; MCEF3_fb1[82] = 0.00638945;//eta is 0.65
3730 MCEF1_fb1[83] = 0.00484419; MCEF3_fb1[83] = 0.0067127;//eta is 0.67
3731 MCEF1_fb1[84] = 0.00484068; MCEF3_fb1[84] = 0.00621438;//eta is 0.69
3732 MCEF1_fb1[85] = 0.00480704; MCEF3_fb1[85] = 0.00606007;//eta is 0.71
3733 MCEF1_fb1[86] = 0.00478109; MCEF3_fb1[86] = 0.00596953;//eta is 0.73
3734 MCEF1_fb1[87] = 0.00488645; MCEF3_fb1[87] = 0.00594765;//eta is 0.75
3735 MCEF1_fb1[88] = 0.00459445; MCEF3_fb1[88] = 0.00579464;//eta is 0.77
3736 MCEF1_fb1[89] = 0.00463147; MCEF3_fb1[89] = 0.0056488;//eta is 0.79
3737 MCEF1_fb1[90] = 0.00464974; MCEF3_fb1[90] = 0.00567165;//eta is 0.81
3738 MCEF1_fb1[91] = 0.00458123; MCEF3_fb1[91] = 0.00532498;//eta is 0.83
3739 MCEF1_fb1[92] = 0.0047586; MCEF3_fb1[92] = 0.00472428;//eta is 0.85
3740 MCEF1_fb1[93] = 0.00492789; MCEF3_fb1[93] = 0.00514747;//eta is 0.87
3741 MCEF1_fb1[94] = 0.00464254; MCEF3_fb1[94] = 0.00544923;//eta is 0.89
3742 MCEF1_fb1[95] = 0.00453585; MCEF3_fb1[95] = 0.00511096;//eta is 0.91
3743 MCEF1_fb1[96] = 0.00427511; MCEF3_fb1[96] = 0.00556982;//eta is 0.93
3744 MCEF1_fb1[97] = 0.00463089; MCEF3_fb1[97] = 0.00584684;//eta is 0.95
3745 MCEF1_fb1[98] = 0.00446822; MCEF3_fb1[98] = 0.00615053;//eta is 0.97
3746 MCEF1_fb1[99] = 0.00433282; MCEF3_fb1[99] = 0.0063621;//eta is 0.99
3747
3748 // --- filter bit 128 (TPC only tracks)
3749 MCEF1_fb128[0] = 0.00591299; MCEF3_fb128[0] = 0.00383748;//eta is -0.99
3750 MCEF1_fb128[1] = 0.00588992; MCEF3_fb128[1] = 0.003667;//eta is -0.97
3751 MCEF1_fb128[2] = 0.00570083; MCEF3_fb128[2] = 0.00361253;//eta is -0.95
3752 MCEF1_fb128[3] = 0.00560151; MCEF3_fb128[3] = 0.00333053;//eta is -0.93
3753 MCEF1_fb128[4] = 0.00576251; MCEF3_fb128[4] = 0.00340385;//eta is -0.91
3754 MCEF1_fb128[5] = 0.00554351; MCEF3_fb128[5] = 0.00292592;//eta is -0.89
3755 MCEF1_fb128[6] = 0.0056127; MCEF3_fb128[6] = 0.00332485;//eta is -0.87
3756 MCEF1_fb128[7] = 0.00541225; MCEF3_fb128[7] = 0.00331356;//eta is -0.85
3757 MCEF1_fb128[8] = 0.00525774; MCEF3_fb128[8] = 0.00370341;//eta is -0.83
3758 MCEF1_fb128[9] = 0.00535911; MCEF3_fb128[9] = 0.00380439;//eta is -0.81
3759 MCEF1_fb128[10] = 0.00532029; MCEF3_fb128[10] = 0.00418114;//eta is -0.79
3760 MCEF1_fb128[11] = 0.00542025; MCEF3_fb128[11] = 0.00387263;//eta is -0.77
3761 MCEF1_fb128[12] = 0.00520206; MCEF3_fb128[12] = 0.00414251;//eta is -0.75
3762 MCEF1_fb128[13] = 0.00513343; MCEF3_fb128[13] = 0.00400127;//eta is -0.73
3763 MCEF1_fb128[14] = 0.00543818; MCEF3_fb128[14] = 0.00416512;//eta is -0.71
3764 MCEF1_fb128[15] = 0.00553948; MCEF3_fb128[15] = 0.00406533;//eta is -0.69
3765 MCEF1_fb128[16] = 0.00556544; MCEF3_fb128[16] = 0.00428908;//eta is -0.67
3766 MCEF1_fb128[17] = 0.00558886; MCEF3_fb128[17] = 0.00435018;//eta is -0.65
3767 MCEF1_fb128[18] = 0.0057369; MCEF3_fb128[18] = 0.00445061;//eta is -0.63
3768 MCEF1_fb128[19] = 0.00616797; MCEF3_fb128[19] = 0.00443836;//eta is -0.61
3769 MCEF1_fb128[20] = 0.00604426; MCEF3_fb128[20] = 0.00466612;//eta is -0.59
3770 MCEF1_fb128[21] = 0.00581836; MCEF3_fb128[21] = 0.00473329;//eta is -0.57
3771 MCEF1_fb128[22] = 0.00573374; MCEF3_fb128[22] = 0.00429524;//eta is -0.55
3772 MCEF1_fb128[23] = 0.00563837; MCEF3_fb128[23] = 0.0047092;//eta is -0.53
3773 MCEF1_fb128[24] = 0.00628575; MCEF3_fb128[24] = 0.0050631;//eta is -0.51
3774 MCEF1_fb128[25] = 0.00606812; MCEF3_fb128[25] = 0.00480595;//eta is -0.49
3775 MCEF1_fb128[26] = 0.00614111; MCEF3_fb128[26] = 0.00474149;//eta is -0.47
3776 MCEF1_fb128[27] = 0.0060786; MCEF3_fb128[27] = 0.00478341;//eta is -0.45
3777 MCEF1_fb128[28] = 0.00621897; MCEF3_fb128[28] = 0.00499422;//eta is -0.43
3778 MCEF1_fb128[29] = 0.00633416; MCEF3_fb128[29] = 0.00498704;//eta is -0.41
3779 MCEF1_fb128[30] = 0.00638646; MCEF3_fb128[30] = 0.00503349;//eta is -0.39
3780 MCEF1_fb128[31] = 0.00658823; MCEF3_fb128[31] = 0.00485434;//eta is -0.37
3781 MCEF1_fb128[32] = 0.0064771; MCEF3_fb128[32] = 0.00482841;//eta is -0.35
3782 MCEF1_fb128[33] = 0.00694801; MCEF3_fb128[33] = 0.00471101;//eta is -0.33
3783 MCEF1_fb128[34] = 0.00636795; MCEF3_fb128[34] = 0.00491829;//eta is -0.31
3784 MCEF1_fb128[35] = 0.00650891; MCEF3_fb128[35] = 0.00498886;//eta is -0.29
3785 MCEF1_fb128[36] = 0.00653225; MCEF3_fb128[36] = 0.00481843;//eta is -0.27
3786 MCEF1_fb128[37] = 0.00657161; MCEF3_fb128[37] = 0.00511179;//eta is -0.25
3787 MCEF1_fb128[38] = 0.00640258; MCEF3_fb128[38] = 0.00511465;//eta is -0.23
3788 MCEF1_fb128[39] = 0.00626697; MCEF3_fb128[39] = 0.00563881;//eta is -0.21
3789 MCEF1_fb128[40] = 0.00691174; MCEF3_fb128[40] = 0.00534971;//eta is -0.19
3790 MCEF1_fb128[41] = 0.00620562; MCEF3_fb128[41] = 0.00557349;//eta is -0.17
3791 MCEF1_fb128[42] = 0.00649612; MCEF3_fb128[42] = 0.00498584;//eta is -0.15
3792 MCEF1_fb128[43] = 0.0065798; MCEF3_fb128[43] = 0.00542058;//eta is -0.13
3793 MCEF1_fb128[44] = 0.00642303; MCEF3_fb128[44] = 0.00613666;//eta is -0.11
3794 MCEF1_fb128[45] = 0.00653218; MCEF3_fb128[45] = 0.00546797;//eta is -0.09
3795 MCEF1_fb128[46] = 0.0064394; MCEF3_fb128[46] = 0.00555092;//eta is -0.07
3796 MCEF1_fb128[47] = 0.00622478; MCEF3_fb128[47] = 0.00629308;//eta is -0.05
3797 MCEF1_fb128[48] = 0.00639885; MCEF3_fb128[48] = 0.0063706;//eta is -0.03
3798 MCEF1_fb128[49] = 0.00666821; MCEF3_fb128[49] = 0.00631535;//eta is -0.01
3799 MCEF1_fb128[50] = 0.00653648; MCEF3_fb128[50] = 0.00710942;//eta is 0.01
3800 MCEF1_fb128[51] = 0.00638104; MCEF3_fb128[51] = 0.00709238;//eta is 0.03
3801 MCEF1_fb128[52] = 0.0066302; MCEF3_fb128[52] = 0.00752929;//eta is 0.05
3802 MCEF1_fb128[53] = 0.00660631; MCEF3_fb128[53] = 0.00756394;//eta is 0.07
3803 MCEF1_fb128[54] = 0.0068361; MCEF3_fb128[54] = 0.00763778;//eta is 0.09
3804 MCEF1_fb128[55] = 0.00680126; MCEF3_fb128[55] = 0.00766098;//eta is 0.11
3805 MCEF1_fb128[56] = 0.00673586; MCEF3_fb128[56] = 0.0077562;//eta is 0.13
3806 MCEF1_fb128[57] = 0.0066831; MCEF3_fb128[57] = 0.00749726;//eta is 0.15
3807 MCEF1_fb128[58] = 0.00632194; MCEF3_fb128[58] = 0.00814879;//eta is 0.17
3808 MCEF1_fb128[59] = 0.00615406; MCEF3_fb128[59] = 0.00773741;//eta is 0.19
3809 MCEF1_fb128[60] = 0.00640291; MCEF3_fb128[60] = 0.00709483;//eta is 0.21
3810 MCEF1_fb128[61] = 0.00626011; MCEF3_fb128[61] = 0.00748346;//eta is 0.23
3811 MCEF1_fb128[62] = 0.00588407; MCEF3_fb128[62] = 0.00746276;//eta is 0.25
3812 MCEF1_fb128[63] = 0.00615514; MCEF3_fb128[63] = 0.0073329;//eta is 0.27
3813 MCEF1_fb128[64] = 0.00581758; MCEF3_fb128[64] = 0.00749753;//eta is 0.29
3814 MCEF1_fb128[65] = 0.00552913; MCEF3_fb128[65] = 0.00758203;//eta is 0.31
3815 MCEF1_fb128[66] = 0.00618476; MCEF3_fb128[66] = 0.00736015;//eta is 0.33
3816 MCEF1_fb128[67] = 0.00571961; MCEF3_fb128[67] = 0.00717504;//eta is 0.35
3817 MCEF1_fb128[68] = 0.00596315; MCEF3_fb128[68] = 0.00731385;//eta is 0.37
3818 MCEF1_fb128[69] = 0.00617447; MCEF3_fb128[69] = 0.00715388;//eta is 0.39
3819 MCEF1_fb128[70] = 0.00548738; MCEF3_fb128[70] = 0.00691641;//eta is 0.41
3820 MCEF1_fb128[71] = 0.00526914; MCEF3_fb128[71] = 0.00712093;//eta is 0.43
3821 MCEF1_fb128[72] = 0.00549674; MCEF3_fb128[72] = 0.00689862;//eta is 0.45
3822 MCEF1_fb128[73] = 0.00556518; MCEF3_fb128[73] = 0.00701435;//eta is 0.47
3823 MCEF1_fb128[74] = 0.005665; MCEF3_fb128[74] = 0.00676622;//eta is 0.49
3824 MCEF1_fb128[75] = 0.00552991; MCEF3_fb128[75] = 0.00679297;//eta is 0.51
3825 MCEF1_fb128[76] = 0.00569837; MCEF3_fb128[76] = 0.00662394;//eta is 0.53
3826 MCEF1_fb128[77] = 0.00562372; MCEF3_fb128[77] = 0.00638628;//eta is 0.55
3827 MCEF1_fb128[78] = 0.00543001; MCEF3_fb128[78] = 0.00657229;//eta is 0.57
3828 MCEF1_fb128[79] = 0.00524306; MCEF3_fb128[79] = 0.00635786;//eta is 0.59
3829 MCEF1_fb128[80] = 0.00535293; MCEF3_fb128[80] = 0.00641706;//eta is 0.61
3830 MCEF1_fb128[81] = 0.00514873; MCEF3_fb128[81] = 0.00636306;//eta is 0.63
3831 MCEF1_fb128[82] = 0.00510695; MCEF3_fb128[82] = 0.00638945;//eta is 0.65
3832 MCEF1_fb128[83] = 0.00484419; MCEF3_fb128[83] = 0.0067127;//eta is 0.67
3833 MCEF1_fb128[84] = 0.00484068; MCEF3_fb128[84] = 0.00621438;//eta is 0.69
3834 MCEF1_fb128[85] = 0.00480704; MCEF3_fb128[85] = 0.00606007;//eta is 0.71
3835 MCEF1_fb128[86] = 0.00478109; MCEF3_fb128[86] = 0.00596953;//eta is 0.73
3836 MCEF1_fb128[87] = 0.00488645; MCEF3_fb128[87] = 0.00594765;//eta is 0.75
3837 MCEF1_fb128[88] = 0.00459445; MCEF3_fb128[88] = 0.00579464;//eta is 0.77
3838 MCEF1_fb128[89] = 0.00463147; MCEF3_fb128[89] = 0.0056488;//eta is 0.79
3839 MCEF1_fb128[90] = 0.00464974; MCEF3_fb128[90] = 0.00567165;//eta is 0.81
3840 MCEF1_fb128[91] = 0.00458123; MCEF3_fb128[91] = 0.00532498;//eta is 0.83
3841 MCEF1_fb128[92] = 0.0047586; MCEF3_fb128[92] = 0.00472428;//eta is 0.85
3842 MCEF1_fb128[93] = 0.00492789; MCEF3_fb128[93] = 0.00514747;//eta is 0.87
3843 MCEF1_fb128[94] = 0.00464254; MCEF3_fb128[94] = 0.00544923;//eta is 0.89
3844 MCEF1_fb128[95] = 0.00453585; MCEF3_fb128[95] = 0.00511096;//eta is 0.91
3845 MCEF1_fb128[96] = 0.00427511; MCEF3_fb128[96] = 0.00556982;//eta is 0.93
3846 MCEF1_fb128[97] = 0.00463089; MCEF3_fb128[97] = 0.00584684;//eta is 0.95
3847 MCEF1_fb128[98] = 0.00446822; MCEF3_fb128[98] = 0.00615053;//eta is 0.97
3848 MCEF1_fb128[99] = 0.00433282; MCEF3_fb128[99] = 0.0063621;//eta is 0.99
3849
3850 // --- filter bit 272 (hybrid tracks)
3851 MCEF1_fb272[0] = 0.00591299; MCEF3_fb272[0] = 0.00383748;//eta is -0.99
3852 MCEF1_fb272[1] = 0.00588992; MCEF3_fb272[1] = 0.003667;//eta is -0.97
3853 MCEF1_fb272[2] = 0.00570083; MCEF3_fb272[2] = 0.00361253;//eta is -0.95
3854 MCEF1_fb272[3] = 0.00560151; MCEF3_fb272[3] = 0.00333053;//eta is -0.93
3855 MCEF1_fb272[4] = 0.00576251; MCEF3_fb272[4] = 0.00340385;//eta is -0.91
3856 MCEF1_fb272[5] = 0.00554351; MCEF3_fb272[5] = 0.00292592;//eta is -0.89
3857 MCEF1_fb272[6] = 0.0056127; MCEF3_fb272[6] = 0.00332485;//eta is -0.87
3858 MCEF1_fb272[7] = 0.00541225; MCEF3_fb272[7] = 0.00331356;//eta is -0.85
3859 MCEF1_fb272[8] = 0.00525774; MCEF3_fb272[8] = 0.00370341;//eta is -0.83
3860 MCEF1_fb272[9] = 0.00535911; MCEF3_fb272[9] = 0.00380439;//eta is -0.81
3861 MCEF1_fb272[10] = 0.00532029; MCEF3_fb272[10] = 0.00418114;//eta is -0.79
3862 MCEF1_fb272[11] = 0.00542025; MCEF3_fb272[11] = 0.00387263;//eta is -0.77
3863 MCEF1_fb272[12] = 0.00520206; MCEF3_fb272[12] = 0.00414251;//eta is -0.75
3864 MCEF1_fb272[13] = 0.00513343; MCEF3_fb272[13] = 0.00400127;//eta is -0.73
3865 MCEF1_fb272[14] = 0.00543818; MCEF3_fb272[14] = 0.00416512;//eta is -0.71
3866 MCEF1_fb272[15] = 0.00553948; MCEF3_fb272[15] = 0.00406533;//eta is -0.69
3867 MCEF1_fb272[16] = 0.00556544; MCEF3_fb272[16] = 0.00428908;//eta is -0.67
3868 MCEF1_fb272[17] = 0.00558886; MCEF3_fb272[17] = 0.00435018;//eta is -0.65
3869 MCEF1_fb272[18] = 0.0057369; MCEF3_fb272[18] = 0.00445061;//eta is -0.63
3870 MCEF1_fb272[19] = 0.00616797; MCEF3_fb272[19] = 0.00443836;//eta is -0.61
3871 MCEF1_fb272[20] = 0.00604426; MCEF3_fb272[20] = 0.00466612;//eta is -0.59
3872 MCEF1_fb272[21] = 0.00581836; MCEF3_fb272[21] = 0.00473329;//eta is -0.57
3873 MCEF1_fb272[22] = 0.00573374; MCEF3_fb272[22] = 0.00429524;//eta is -0.55
3874 MCEF1_fb272[23] = 0.00563837; MCEF3_fb272[23] = 0.0047092;//eta is -0.53
3875 MCEF1_fb272[24] = 0.00628575; MCEF3_fb272[24] = 0.0050631;//eta is -0.51
3876 MCEF1_fb272[25] = 0.00606812; MCEF3_fb272[25] = 0.00480595;//eta is -0.49
3877 MCEF1_fb272[26] = 0.00614111; MCEF3_fb272[26] = 0.00474149;//eta is -0.47
3878 MCEF1_fb272[27] = 0.0060786; MCEF3_fb272[27] = 0.00478341;//eta is -0.45
3879 MCEF1_fb272[28] = 0.00621897; MCEF3_fb272[28] = 0.00499422;//eta is -0.43
3880 MCEF1_fb272[29] = 0.00633416; MCEF3_fb272[29] = 0.00498704;//eta is -0.41
3881 MCEF1_fb272[30] = 0.00638646; MCEF3_fb272[30] = 0.00503349;//eta is -0.39
3882 MCEF1_fb272[31] = 0.00658823; MCEF3_fb272[31] = 0.00485434;//eta is -0.37
3883 MCEF1_fb272[32] = 0.0064771; MCEF3_fb272[32] = 0.00482841;//eta is -0.35
3884 MCEF1_fb272[33] = 0.00694801; MCEF3_fb272[33] = 0.00471101;//eta is -0.33
3885 MCEF1_fb272[34] = 0.00636795; MCEF3_fb272[34] = 0.00491829;//eta is -0.31
3886 MCEF1_fb272[35] = 0.00650891; MCEF3_fb272[35] = 0.00498886;//eta is -0.29
3887 MCEF1_fb272[36] = 0.00653225; MCEF3_fb272[36] = 0.00481843;//eta is -0.27
3888 MCEF1_fb272[37] = 0.00657161; MCEF3_fb272[37] = 0.00511179;//eta is -0.25
3889 MCEF1_fb272[38] = 0.00640258; MCEF3_fb272[38] = 0.00511465;//eta is -0.23
3890 MCEF1_fb272[39] = 0.00626697; MCEF3_fb272[39] = 0.00563881;//eta is -0.21
3891 MCEF1_fb272[40] = 0.00691174; MCEF3_fb272[40] = 0.00534971;//eta is -0.19
3892 MCEF1_fb272[41] = 0.00620562; MCEF3_fb272[41] = 0.00557349;//eta is -0.17
3893 MCEF1_fb272[42] = 0.00649612; MCEF3_fb272[42] = 0.00498584;//eta is -0.15
3894 MCEF1_fb272[43] = 0.0065798; MCEF3_fb272[43] = 0.00542058;//eta is -0.13
3895 MCEF1_fb272[44] = 0.00642303; MCEF3_fb272[44] = 0.00613666;//eta is -0.11
3896 MCEF1_fb272[45] = 0.00653218; MCEF3_fb272[45] = 0.00546797;//eta is -0.09
3897 MCEF1_fb272[46] = 0.0064394; MCEF3_fb272[46] = 0.00555092;//eta is -0.07
3898 MCEF1_fb272[47] = 0.00622478; MCEF3_fb272[47] = 0.00629308;//eta is -0.05
3899 MCEF1_fb272[48] = 0.00639885; MCEF3_fb272[48] = 0.0063706;//eta is -0.03
3900 MCEF1_fb272[49] = 0.00666821; MCEF3_fb272[49] = 0.00631535;//eta is -0.01
3901 MCEF1_fb272[50] = 0.00653648; MCEF3_fb272[50] = 0.00710942;//eta is 0.01
3902 MCEF1_fb272[51] = 0.00638104; MCEF3_fb272[51] = 0.00709238;//eta is 0.03
3903 MCEF1_fb272[52] = 0.0066302; MCEF3_fb272[52] = 0.00752929;//eta is 0.05
3904 MCEF1_fb272[53] = 0.00660631; MCEF3_fb272[53] = 0.00756394;//eta is 0.07
3905 MCEF1_fb272[54] = 0.0068361; MCEF3_fb272[54] = 0.00763778;//eta is 0.09
3906 MCEF1_fb272[55] = 0.00680126; MCEF3_fb272[55] = 0.00766098;//eta is 0.11
3907 MCEF1_fb272[56] = 0.00673586; MCEF3_fb272[56] = 0.0077562;//eta is 0.13
3908 MCEF1_fb272[57] = 0.0066831; MCEF3_fb272[57] = 0.00749726;//eta is 0.15
3909 MCEF1_fb272[58] = 0.00632194; MCEF3_fb272[58] = 0.00814879;//eta is 0.17
3910 MCEF1_fb272[59] = 0.00615406; MCEF3_fb272[59] = 0.00773741;//eta is 0.19
3911 MCEF1_fb272[60] = 0.00640291; MCEF3_fb272[60] = 0.00709483;//eta is 0.21
3912 MCEF1_fb272[61] = 0.00626011; MCEF3_fb272[61] = 0.00748346;//eta is 0.23
3913 MCEF1_fb272[62] = 0.00588407; MCEF3_fb272[62] = 0.00746276;//eta is 0.25
3914 MCEF1_fb272[63] = 0.00615514; MCEF3_fb272[63] = 0.0073329;//eta is 0.27
3915 MCEF1_fb272[64] = 0.00581758; MCEF3_fb272[64] = 0.00749753;//eta is 0.29
3916 MCEF1_fb272[65] = 0.00552913; MCEF3_fb272[65] = 0.00758203;//eta is 0.31
3917 MCEF1_fb272[66] = 0.00618476; MCEF3_fb272[66] = 0.00736015;//eta is 0.33
3918 MCEF1_fb272[67] = 0.00571961; MCEF3_fb272[67] = 0.00717504;//eta is 0.35
3919 MCEF1_fb272[68] = 0.00596315; MCEF3_fb272[68] = 0.00731385;//eta is 0.37
3920 MCEF1_fb272[69] = 0.00617447; MCEF3_fb272[69] = 0.00715388;//eta is 0.39
3921 MCEF1_fb272[70] = 0.00548738; MCEF3_fb272[70] = 0.00691641;//eta is 0.41
3922 MCEF1_fb272[71] = 0.00526914; MCEF3_fb272[71] = 0.00712093;//eta is 0.43
3923 MCEF1_fb272[72] = 0.00549674; MCEF3_fb272[72] = 0.00689862;//eta is 0.45
3924 MCEF1_fb272[73] = 0.00556518; MCEF3_fb272[73] = 0.00701435;//eta is 0.47
3925 MCEF1_fb272[74] = 0.005665; MCEF3_fb272[74] = 0.00676622;//eta is 0.49
3926 MCEF1_fb272[75] = 0.00552991; MCEF3_fb272[75] = 0.00679297;//eta is 0.51
3927 MCEF1_fb272[76] = 0.00569837; MCEF3_fb272[76] = 0.00662394;//eta is 0.53
3928 MCEF1_fb272[77] = 0.00562372; MCEF3_fb272[77] = 0.00638628;//eta is 0.55
3929 MCEF1_fb272[78] = 0.00543001; MCEF3_fb272[78] = 0.00657229;//eta is 0.57
3930 MCEF1_fb272[79] = 0.00524306; MCEF3_fb272[79] = 0.00635786;//eta is 0.59
3931 MCEF1_fb272[80] = 0.00535293; MCEF3_fb272[80] = 0.00641706;//eta is 0.61
3932 MCEF1_fb272[81] = 0.00514873; MCEF3_fb272[81] = 0.00636306;//eta is 0.63
3933 MCEF1_fb272[82] = 0.00510695; MCEF3_fb272[82] = 0.00638945;//eta is 0.65
3934 MCEF1_fb272[83] = 0.00484419; MCEF3_fb272[83] = 0.0067127;//eta is 0.67
3935 MCEF1_fb272[84] = 0.00484068; MCEF3_fb272[84] = 0.00621438;//eta is 0.69
3936 MCEF1_fb272[85] = 0.00480704; MCEF3_fb272[85] = 0.00606007;//eta is 0.71
3937 MCEF1_fb272[86] = 0.00478109; MCEF3_fb272[86] = 0.00596953;//eta is 0.73
3938 MCEF1_fb272[87] = 0.00488645; MCEF3_fb272[87] = 0.00594765;//eta is 0.75
3939 MCEF1_fb272[88] = 0.00459445; MCEF3_fb272[88] = 0.00579464;//eta is 0.77
3940 MCEF1_fb272[89] = 0.00463147; MCEF3_fb272[89] = 0.0056488;//eta is 0.79
3941 MCEF1_fb272[90] = 0.00464974; MCEF3_fb272[90] = 0.00567165;//eta is 0.81
3942 MCEF1_fb272[91] = 0.00458123; MCEF3_fb272[91] = 0.00532498;//eta is 0.83
3943 MCEF1_fb272[92] = 0.0047586; MCEF3_fb272[92] = 0.00472428;//eta is 0.85
3944 MCEF1_fb272[93] = 0.00492789; MCEF3_fb272[93] = 0.00514747;//eta is 0.87
3945 MCEF1_fb272[94] = 0.00464254; MCEF3_fb272[94] = 0.00544923;//eta is 0.89
3946 MCEF1_fb272[95] = 0.00453585; MCEF3_fb272[95] = 0.00511096;//eta is 0.91
3947 MCEF1_fb272[96] = 0.00427511; MCEF3_fb272[96] = 0.00556982;//eta is 0.93
3948 MCEF1_fb272[97] = 0.00463089; MCEF3_fb272[97] = 0.00584684;//eta is 0.95
3949 MCEF1_fb272[98] = 0.00446822; MCEF3_fb272[98] = 0.00615053;//eta is 0.97
3950 MCEF1_fb272[99] = 0.00433282; MCEF3_fb272[99] = 0.0063621;//eta is 0.99
3951
3952 // --- filter bit 1 with tight DCA cut...
3953 MCEF1_dca[0] = 0.00591299; MCEF3_dca[0] = 0.00383748;//eta is -0.99
3954 MCEF1_dca[1] = 0.00588992; MCEF3_dca[1] = 0.003667;//eta is -0.97
3955 MCEF1_dca[2] = 0.00570083; MCEF3_dca[2] = 0.00361253;//eta is -0.95
3956 MCEF1_dca[3] = 0.00560151; MCEF3_dca[3] = 0.00333053;//eta is -0.93
3957 MCEF1_dca[4] = 0.00576251; MCEF3_dca[4] = 0.00340385;//eta is -0.91
3958 MCEF1_dca[5] = 0.00554351; MCEF3_dca[5] = 0.00292592;//eta is -0.89
3959 MCEF1_dca[6] = 0.0056127; MCEF3_dca[6] = 0.00332485;//eta is -0.87
3960 MCEF1_dca[7] = 0.00541225; MCEF3_dca[7] = 0.00331356;//eta is -0.85
3961 MCEF1_dca[8] = 0.00525774; MCEF3_dca[8] = 0.00370341;//eta is -0.83
3962 MCEF1_dca[9] = 0.00535911; MCEF3_dca[9] = 0.00380439;//eta is -0.81
3963 MCEF1_dca[10] = 0.00532029; MCEF3_dca[10] = 0.00418114;//eta is -0.79
3964 MCEF1_dca[11] = 0.00542025; MCEF3_dca[11] = 0.00387263;//eta is -0.77
3965 MCEF1_dca[12] = 0.00520206; MCEF3_dca[12] = 0.00414251;//eta is -0.75
3966 MCEF1_dca[13] = 0.00513343; MCEF3_dca[13] = 0.00400127;//eta is -0.73
3967 MCEF1_dca[14] = 0.00543818; MCEF3_dca[14] = 0.00416512;//eta is -0.71
3968 MCEF1_dca[15] = 0.00553948; MCEF3_dca[15] = 0.00406533;//eta is -0.69
3969 MCEF1_dca[16] = 0.00556544; MCEF3_dca[16] = 0.00428908;//eta is -0.67
3970 MCEF1_dca[17] = 0.00558886; MCEF3_dca[17] = 0.00435018;//eta is -0.65
3971 MCEF1_dca[18] = 0.0057369; MCEF3_dca[18] = 0.00445061;//eta is -0.63
3972 MCEF1_dca[19] = 0.00616797; MCEF3_dca[19] = 0.00443836;//eta is -0.61
3973 MCEF1_dca[20] = 0.00604426; MCEF3_dca[20] = 0.00466612;//eta is -0.59
3974 MCEF1_dca[21] = 0.00581836; MCEF3_dca[21] = 0.00473329;//eta is -0.57
3975 MCEF1_dca[22] = 0.00573374; MCEF3_dca[22] = 0.00429524;//eta is -0.55
3976 MCEF1_dca[23] = 0.00563837; MCEF3_dca[23] = 0.0047092;//eta is -0.53
3977 MCEF1_dca[24] = 0.00628575; MCEF3_dca[24] = 0.0050631;//eta is -0.51
3978 MCEF1_dca[25] = 0.00606812; MCEF3_dca[25] = 0.00480595;//eta is -0.49
3979 MCEF1_dca[26] = 0.00614111; MCEF3_dca[26] = 0.00474149;//eta is -0.47
3980 MCEF1_dca[27] = 0.0060786; MCEF3_dca[27] = 0.00478341;//eta is -0.45
3981 MCEF1_dca[28] = 0.00621897; MCEF3_dca[28] = 0.00499422;//eta is -0.43
3982 MCEF1_dca[29] = 0.00633416; MCEF3_dca[29] = 0.00498704;//eta is -0.41
3983 MCEF1_dca[30] = 0.00638646; MCEF3_dca[30] = 0.00503349;//eta is -0.39
3984 MCEF1_dca[31] = 0.00658823; MCEF3_dca[31] = 0.00485434;//eta is -0.37
3985 MCEF1_dca[32] = 0.0064771; MCEF3_dca[32] = 0.00482841;//eta is -0.35
3986 MCEF1_dca[33] = 0.00694801; MCEF3_dca[33] = 0.00471101;//eta is -0.33
3987 MCEF1_dca[34] = 0.00636795; MCEF3_dca[34] = 0.00491829;//eta is -0.31
3988 MCEF1_dca[35] = 0.00650891; MCEF3_dca[35] = 0.00498886;//eta is -0.29
3989 MCEF1_dca[36] = 0.00653225; MCEF3_dca[36] = 0.00481843;//eta is -0.27
3990 MCEF1_dca[37] = 0.00657161; MCEF3_dca[37] = 0.00511179;//eta is -0.25
3991 MCEF1_dca[38] = 0.00640258; MCEF3_dca[38] = 0.00511465;//eta is -0.23
3992 MCEF1_dca[39] = 0.00626697; MCEF3_dca[39] = 0.00563881;//eta is -0.21
3993 MCEF1_dca[40] = 0.00691174; MCEF3_dca[40] = 0.00534971;//eta is -0.19
3994 MCEF1_dca[41] = 0.00620562; MCEF3_dca[41] = 0.00557349;//eta is -0.17
3995 MCEF1_dca[42] = 0.00649612; MCEF3_dca[42] = 0.00498584;//eta is -0.15
3996 MCEF1_dca[43] = 0.0065798; MCEF3_dca[43] = 0.00542058;//eta is -0.13
3997 MCEF1_dca[44] = 0.00642303; MCEF3_dca[44] = 0.00613666;//eta is -0.11
3998 MCEF1_dca[45] = 0.00653218; MCEF3_dca[45] = 0.00546797;//eta is -0.09
3999 MCEF1_dca[46] = 0.0064394; MCEF3_dca[46] = 0.00555092;//eta is -0.07
4000 MCEF1_dca[47] = 0.00622478; MCEF3_dca[47] = 0.00629308;//eta is -0.05
4001 MCEF1_dca[48] = 0.00639885; MCEF3_dca[48] = 0.0063706;//eta is -0.03
4002 MCEF1_dca[49] = 0.00666821; MCEF3_dca[49] = 0.00631535;//eta is -0.01
4003 MCEF1_dca[50] = 0.00653648; MCEF3_dca[50] = 0.00710942;//eta is 0.01
4004 MCEF1_dca[51] = 0.00638104; MCEF3_dca[51] = 0.00709238;//eta is 0.03
4005 MCEF1_dca[52] = 0.0066302; MCEF3_dca[52] = 0.00752929;//eta is 0.05
4006 MCEF1_dca[53] = 0.00660631; MCEF3_dca[53] = 0.00756394;//eta is 0.07
4007 MCEF1_dca[54] = 0.0068361; MCEF3_dca[54] = 0.00763778;//eta is 0.09
4008 MCEF1_dca[55] = 0.00680126; MCEF3_dca[55] = 0.00766098;//eta is 0.11
4009 MCEF1_dca[56] = 0.00673586; MCEF3_dca[56] = 0.0077562;//eta is 0.13
4010 MCEF1_dca[57] = 0.0066831; MCEF3_dca[57] = 0.00749726;//eta is 0.15
4011 MCEF1_dca[58] = 0.00632194; MCEF3_dca[58] = 0.00814879;//eta is 0.17
4012 MCEF1_dca[59] = 0.00615406; MCEF3_dca[59] = 0.00773741;//eta is 0.19
4013 MCEF1_dca[60] = 0.00640291; MCEF3_dca[60] = 0.00709483;//eta is 0.21
4014 MCEF1_dca[61] = 0.00626011; MCEF3_dca[61] = 0.00748346;//eta is 0.23
4015 MCEF1_dca[62] = 0.00588407; MCEF3_dca[62] = 0.00746276;//eta is 0.25
4016 MCEF1_dca[63] = 0.00615514; MCEF3_dca[63] = 0.0073329;//eta is 0.27
4017 MCEF1_dca[64] = 0.00581758; MCEF3_dca[64] = 0.00749753;//eta is 0.29
4018 MCEF1_dca[65] = 0.00552913; MCEF3_dca[65] = 0.00758203;//eta is 0.31
4019 MCEF1_dca[66] = 0.00618476; MCEF3_dca[66] = 0.00736015;//eta is 0.33
4020 MCEF1_dca[67] = 0.00571961; MCEF3_dca[67] = 0.00717504;//eta is 0.35
4021 MCEF1_dca[68] = 0.00596315; MCEF3_dca[68] = 0.00731385;//eta is 0.37
4022 MCEF1_dca[69] = 0.00617447; MCEF3_dca[69] = 0.00715388;//eta is 0.39
4023 MCEF1_dca[70] = 0.00548738; MCEF3_dca[70] = 0.00691641;//eta is 0.41
4024 MCEF1_dca[71] = 0.00526914; MCEF3_dca[71] = 0.00712093;//eta is 0.43
4025 MCEF1_dca[72] = 0.00549674; MCEF3_dca[72] = 0.00689862;//eta is 0.45
4026 MCEF1_dca[73] = 0.00556518; MCEF3_dca[73] = 0.00701435;//eta is 0.47
4027 MCEF1_dca[74] = 0.005665; MCEF3_dca[74] = 0.00676622;//eta is 0.49
4028 MCEF1_dca[75] = 0.00552991; MCEF3_dca[75] = 0.00679297;//eta is 0.51
4029 MCEF1_dca[76] = 0.00569837; MCEF3_dca[76] = 0.00662394;//eta is 0.53
4030 MCEF1_dca[77] = 0.00562372; MCEF3_dca[77] = 0.00638628;//eta is 0.55
4031 MCEF1_dca[78] = 0.00543001; MCEF3_dca[78] = 0.00657229;//eta is 0.57
4032 MCEF1_dca[79] = 0.00524306; MCEF3_dca[79] = 0.00635786;//eta is 0.59
4033 MCEF1_dca[80] = 0.00535293; MCEF3_dca[80] = 0.00641706;//eta is 0.61
4034 MCEF1_dca[81] = 0.00514873; MCEF3_dca[81] = 0.00636306;//eta is 0.63
4035 MCEF1_dca[82] = 0.00510695; MCEF3_dca[82] = 0.00638945;//eta is 0.65
4036 MCEF1_dca[83] = 0.00484419; MCEF3_dca[83] = 0.0067127;//eta is 0.67
4037 MCEF1_dca[84] = 0.00484068; MCEF3_dca[84] = 0.00621438;//eta is 0.69
4038 MCEF1_dca[85] = 0.00480704; MCEF3_dca[85] = 0.00606007;//eta is 0.71
4039 MCEF1_dca[86] = 0.00478109; MCEF3_dca[86] = 0.00596953;//eta is 0.73
4040 MCEF1_dca[87] = 0.00488645; MCEF3_dca[87] = 0.00594765;//eta is 0.75
4041 MCEF1_dca[88] = 0.00459445; MCEF3_dca[88] = 0.00579464;//eta is 0.77
4042 MCEF1_dca[89] = 0.00463147; MCEF3_dca[89] = 0.0056488;//eta is 0.79
4043 MCEF1_dca[90] = 0.00464974; MCEF3_dca[90] = 0.00567165;//eta is 0.81
4044 MCEF1_dca[91] = 0.00458123; MCEF3_dca[91] = 0.00532498;//eta is 0.83
4045 MCEF1_dca[92] = 0.0047586; MCEF3_dca[92] = 0.00472428;//eta is 0.85
4046 MCEF1_dca[93] = 0.00492789; MCEF3_dca[93] = 0.00514747;//eta is 0.87
4047 MCEF1_dca[94] = 0.00464254; MCEF3_dca[94] = 0.00544923;//eta is 0.89
4048 MCEF1_dca[95] = 0.00453585; MCEF3_dca[95] = 0.00511096;//eta is 0.91
4049 MCEF1_dca[96] = 0.00427511; MCEF3_dca[96] = 0.00556982;//eta is 0.93
4050 MCEF1_dca[97] = 0.00463089; MCEF3_dca[97] = 0.00584684;//eta is 0.95
4051 MCEF1_dca[98] = 0.00446822; MCEF3_dca[98] = 0.00615053;//eta is 0.97
4052 MCEF1_dca[99] = 0.00433282; MCEF3_dca[99] = 0.0063621;//eta is 0.99
4053
4054
4055
4056}
4057
4058
4059
4060float AliAnalysisTaskCMEv2A::calc3(float xn, float yn, float x2n, float y2n, float M)
4061{
4062
4063 float Qn2 = xn*xn+yn*yn;
4064 float Q2n2 = x2n*x2n+y2n*y2n;
4065 float Qn2d = xn*xn-yn*yn;
4066
4067 float W_3 = M*(M-1)*(M-2);
4068
4069 float first = ((x2n*Qn2d)-(2*Qn2))/W_3;
4070 float second = (Q2n2-(2*M))/W_3;
4071
4072 return first-second;
4073
4074}
4075
4076
4077
4078float AliAnalysisTaskCMEv2A::calc4(float xn, float yn, float x2n, float y2n, float M)
4079{
4080
4081 float Qn2 = xn*xn+yn*yn;
4082 float Q2n2 = x2n*x2n+y2n*y2n;
4083 float Qn4 = Qn2*Qn2;
4084 float Qn2d = xn*xn-yn*yn;
4085
4086 float W_4 = M*(M-1)*(M-2)*(M-3);
4087
4088 float first = (Qn4+Q2n2-(2*(x2n*Qn2d)))/W_4;
4089 float second = 2*((2*(M-2)*Qn2)-(M*(M-3)))/W_4;
4090
4091 return first-second;
4092
4093}
4094
4095
4096// function definition
4097float GetDPhiStar(float phi1, float pt1, float charge1, float phi2, float pt2, float charge2, float radius, float bSign)
4098{
4099 //
4100 // calculates dphistar
4101 //
4102 float dphistar = phi1 - phi2 - charge1 * bSign * asin(0.075 * radius / pt1) + charge2 * bSign * asin(0.075 * radius / pt2);
4103
4104 //static const double kPi = pi;
4105 float kPi = pi;
4106
4107 // circularity
4108 // if (dphistar > 2 * kPi)
4109 // dphistar -= 2 * kPi;
4110 // if (dphistar < -2 * kPi)
4111 // dphistar += 2 * kPi;
4112
4113 if (dphistar > kPi)
4114 dphistar = kPi * 2 - dphistar;
4115 if (dphistar < -kPi)
4116 dphistar = -kPi * 2 - dphistar;
4117 if (dphistar > kPi) // might look funny but is needed
4118 dphistar = kPi * 2 - dphistar;
4119
4120 return dphistar;
4121
4122}
4123
4124