6f911f35648eefa423a5a4aa9efc2f7b4edf0dc1
[u/mrichter/AliRoot.git] / PWGCF / Correlations / DPhi / AliAnalysisTask3PCorrelations.cxx
1 /**************************************************************************
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3  *                                                                        *
4  * Author: The ALICE Off-line Project.                                    *
5  * Contributors are mentioned in the code where appropriate.              *
6  *                                                                        *
7  * Permission to use, copy, modify and distribute this software and its   *
8  * documentation strictly for non-commercial purposes is hereby granted   *
9  * without fee, provided that the above copyright notice appears in all   *
10  * copies and that both the copyright notice and this permission notice   *
11  * appear in the supporting documentation. The authors make no claims     *
12  * about the suitability of this software for any purpose. It is          *
13  * provided "as is" without express or implied warranty.                  *
14  **************************************************************************/
15
16 /* $Id:$ */
17
18 #include <TROOT.h>
19 #include <TChain.h>
20 #include <TFile.h>
21 #include <TList.h>
22 #include <TMath.h>
23 #include <TTree.h>
24 #include <TH1F.h>
25 #include <TH2F.h>
26 #include <TH3F.h>
27 #include <TProfile.h>
28 #include <TH1D.h>
29 #include <TH2D.h>
30 #include <TH3D.h>
31 #include <TRandom.h>
32 #include "AliAnalysisManager.h"
33
34 #include "AliAODHandler.h"
35 #include "AliAODInputHandler.h"
36 //#include "AliAODMCParticle.h"
37 #include "AliInputEventHandler.h"
38 #include "AliLog.h"
39 #include "AliESDEvent.h"
40 #include "AliESDInputHandler.h"
41 #include "AliMultiplicity.h"
42 #include "AliCentrality.h"
43 #include "AliAnalysisTask3PCorrelations.h"
44
45 ClassImp(AliAnalysisTask3PCorrelations)
46
47 AliAnalysisTask3PCorrelations::AliAnalysisTask3PCorrelations()
48 : AliAnalysisTaskSE(),
49 fAODEvent(0), 
50 fESDEvent(0),             //! ESD Event 
51 fInputHandler(0),
52 _outputHistoList(0),
53 _twoPi         ( 2.0 * 3.1415927),
54 _eventCount    ( 0), 
55 _debugLevel    ( 0),
56 _singlesOnly   ( 0), 
57 _useWeights    ( 0), 
58 _rejectPileup  ( 1), 
59 _rejectPairConversion ( 0), 
60 _vertexZMin           ( -10), 
61 _vertexZMax           (  10), 
62 _vertexXYMin          ( -10),
63 _vertexXYMax          (  10),
64 _centralityMethod     (  4),
65 _centralityMin        (  0.),
66 _centralityMax        (  0.),
67 _dcaZMin              ( -3), 
68 _dcaZMax              (  3.), 
69 _dcaXYMin             ( -3.), 
70 _dcaXYMax             (  3.),
71 _dedxMin              ( 0),
72 _dedxMax              ( 100000),
73 _nClusterMin          ( 70), 
74 _trackFilterBit       ( 128),
75 _field    ( 1.),
76 _nTracks  ( 0 ),
77 _mult0    ( 0 ),
78 _mult1    ( 0 ),
79 _mult2    ( 0 ),
80 _mult3    ( 0 ),
81 _mult4    ( 0 ),
82 _mult5    ( 0 ),
83 _mult6    ( 0 ),
84 arraySize ( 2000),
85 _id_1(0),       
86 _charge_1(0),    
87 _iPhi_1(0),    
88 _pt_1(0),       
89 _px_1(0),      
90 _py_1(0),      
91 _pz_1(0),      
92 _correction_1(0),   
93 _dedx_1(0),   
94 _id_2(0),      
95 _charge_2(0),    
96 _iPhi_2(0),    
97 _pt_2(0),      
98 _px_2(0),      
99 _py_2(0),      
100 _pz_2(0),      
101 _correction_2(0),   
102 _dedx_2(0),   
103 _id_3(0),      
104 _charge_3(0),    
105 _iPhi_3(0),    
106 _pt_3(0),      
107 _px_3(0),      
108 _py_3(0),      
109 _pz_3(0),      
110 _correction_3(0),   
111 _dedx_3(0),   
112
113 _correctionWeight_1P(0),   
114 _correctionWeight_1M(0),   
115 _correctionWeight_2P(0),   
116 _correctionWeight_2M(0),   
117 _correctionWeight_3P(0),   
118 _correctionWeight_3M(0),   
119 _nBins_M0(500),       _min_M0(0),        _max_M0(10000),          _width_M0(20),
120 _nBins_M1(500),       _min_M1(0),        _max_M1(10000),          _width_M1(20),
121 _nBins_M2(500),       _min_M2(0),        _max_M2(10000),          _width_M2(20),
122 _nBins_M3(500),       _min_M3(0),        _max_M3(10000),          _width_M3(20),
123 _nBins_M4(100),       _min_M4(0),        _max_M4(1),              _width_M4(0.01),
124 _nBins_M5(100),       _min_M5(0),        _max_M5(1),              _width_M5(0.01),
125 _nBins_M6(100),       _min_M6(0),        _max_M6(1),              _width_M6(0.01),
126 _nBins_vsM(0),        _min_vsM(0),       _max_vsM(0),             _width_vsM(0),
127 _nBins_vertexZ(40),   _min_vertexZ(-10), _max_vertexZ(10),        _width_vertexZ(0.5),
128
129 _nBins_pt_1(8),       _min_pt_1(2.),    _max_pt_1(10.0),          _width_pt_1(1.0),
130 _nBins_phi_1(72),     _min_phi_1(0),     _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
131 _nBins_eta_1(1),      _min_eta_1(-0.25), _max_eta_1(0.25),        _width_eta_1(0.5),
132 _nBins_etaPhi_1(0), 
133 _nBins_etaPhiPt_1(0),
134 _nBins_zEtaPhiPt_1(0),
135
136 _nBins_pt_2(5),       _min_pt_2(0.5),     _max_pt_2(1.0),          _width_pt_2(0.1),
137 _nBins_phi_2(72),     _min_phi_2(0),      _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
138 _nBins_eta_2(1),      _min_eta_2(-1),     _max_eta_2(-0.7),        _width_eta_2(0.3),
139 _nBins_etaPhi_2(0), 
140 _nBins_etaPhiPt_2(0),
141 _nBins_zEtaPhiPt_2(0),
142
143 _nBins_pt_3(5),       _min_pt_3(0.5),     _max_pt_3(1.0),          _width_pt_3(0.1),
144 _nBins_phi_3(72),     _min_phi_3(0),      _max_phi_3(2.*3.1415927),_width_phi_3(2.*3.1415927/72),
145 _nBins_eta_3(1),      _min_eta_3(0.7),    _max_eta_3(1.),          _width_eta_3(0.3),
146 _nBins_etaPhi_3(0), 
147 _nBins_etaPhiPt_3(0),
148 _nBins_zEtaPhiPt_3(0),
149
150 _nBins_phi_12(0),
151 _nBins_phi_13(0),
152 _nBins_phi_23(0),
153 _nBins_phi_123(0),
154
155
156 _nBins_etaPhi_12(0),
157 _nBins_etaPhi_13(0),
158 _nBins_etaPhi_23(0),
159 _nBins_etaPhi_123(0),
160
161 __n1_1(0),
162 __n1_2(0),
163 __n1_3(0),
164 __n2_12(0),   
165 __n2_13(0),   
166 __n2_23(0),   
167 __n3_123(0),   
168 __n1Nw_1(0),
169 __n1Nw_2(0),
170 __n1Nw_3(0),
171 __n2Nw_12(0),   
172 __n2Nw_13(0),   
173 __n2Nw_23(0),   
174 __n3Nw_123(0),  
175
176 __n1_1_vsPt(0),
177 __n1_1_vsPhi(0), 
178 __n1_1P_vsZEtaPhiPt(0),
179 __n1_1M_vsZEtaPhiPt(0),
180
181 __n1_2_vsPt(0),
182 __n1_2_vsPhi(0), 
183 __n1_2P_vsZEtaPhiPt(0),
184 __n1_2M_vsZEtaPhiPt(0),
185
186 __n1_3_vsPt(0),
187 __n1_3_vsPhi(0), 
188 __n1_3P_vsZEtaPhiPt(0),
189 __n1_3M_vsZEtaPhiPt(0),
190
191 __n2_12_vsPhi(0),
192 __n2_13_vsPhi(0),
193 __n2_23_vsPhi(0),
194 __n3_123_vsPhi(0),
195
196 _weight_1P      ( 0    ),
197 _weight_1M      ( 0    ),
198 _weight_2P      ( 0    ),
199 _weight_2M      ( 0    ),
200 _weight_3P      ( 0    ),
201 _weight_3M      ( 0    ),
202 _eventAccounting ( 0),
203 _m0 ( 0),
204 _m1 ( 0),
205 _m2 ( 0),
206 _m3 ( 0),
207 _m4 ( 0),
208 _m5 ( 0),
209 _m6 ( 0),
210 _vertexZ        ( 0),
211
212 _n1_1_vsPhi     ( 0),
213 _n1_1_vsM          ( 0), // w/ weight
214 _n1Nw_1_vsM        ( 0), // w/o weight
215
216 _n1_1_vsPt      ( 0),         
217 _n1_1P_vsZVsEtaVsPhiVsPt ( 0),
218 _n1_1M_vsZVsEtaVsPhiVsPt ( 0),
219 _dedxVsP_1         ( 0),
220 _corrDedxVsP_1     ( 0),
221 _betaVsP_1         ( 0),
222
223 _n1_2_vsPhi        ( 0),
224 _n1_2_vsM          ( 0),
225 _n1Nw_2_vsM        ( 0),
226 _n1_2_vsPt         ( 0),       
227 _n1_2P_vsZVsEtaVsPhiVsPt ( 0), 
228 _n1_2M_vsZVsEtaVsPhiVsPt ( 0), 
229 _dedxVsP_2         ( 0),
230 _corrDedxVsP_2     ( 0),
231 _betaVsP_2         ( 0),
232
233 _n1_3_vsPhi        ( 0),
234 _n1_3_vsM          ( 0),
235 _n1Nw_3_vsM        ( 0),
236 _n1_3_vsPt         ( 0),       
237 _n1_3P_vsZVsEtaVsPhiVsPt ( 0), 
238 _n1_3M_vsZVsEtaVsPhiVsPt ( 0), 
239 _dedxVsP_3         ( 0),
240 _corrDedxVsP_3     ( 0),
241 _betaVsP_3         ( 0),
242
243 _n2_12_vsPhi    ( 0),  
244 _n2_13_vsPhi    ( 0),  
245 _n2_23_vsPhi    ( 0),  
246 _n3_123_vsPhi   ( 0),  
247
248 _n2_12_vsM         ( 0),        
249 _n2_13_vsM         ( 0),        
250 _n2_23_vsM         ( 0),        
251 _n3_123_vsM        ( 0),        
252
253 _n2Nw_12_vsM       ( 0),        
254 _n2Nw_13_vsM       ( 0),        
255 _n2Nw_23_vsM       ( 0),        
256 _n3Nw_123_vsM      ( 0),        
257
258 n1Name("NA"),
259 n2Name("NA"),
260 n3Name("NA"),
261 n1NwName("NA"),
262 n2NwName("NA"),
263 n3NwName("NA"),
264 n1n1Name("NA"),
265 n1n1n1Name("NA"),
266 n2n1Name("NA"),
267 r1Name("NA"),
268 r2Name("NA"),
269 r3Name("NA"),
270 r2r1Name("NA"),
271 c2Name("NA"),
272 c3Name("NA"),
273 d3Name("NA"),
274 p3Name("NA"),
275
276 intR2Name("NA"),
277 binCorrName("NA"),
278 intBinCorrName("NA"),
279
280 countsName("NA"),
281 part_1_Name("NA"),
282 part_1P_Name("NA"),
283 part_1M_Name("NA"),
284 part_2_Name("NA"),
285 part_2P_Name("NA"),
286 part_2M_Name("NA"),
287 part_3_Name("NA"),
288 part_3P_Name("NA"),
289 part_3M_Name("NA"),
290 pair_12_Name("NA"),
291 pair_13_Name("NA"),
292 pair_23_Name("NA"),
293 triplet_Name("NA"),
294
295 _title_counts("NA"),
296 _title_m0("NA"),
297 _title_m1("NA"),
298 _title_m2("NA"),
299 _title_m3("NA"),
300 _title_m4("NA"),
301 _title_m5("NA"),
302 _title_m6("NA"),
303
304 _title_phi_1("NA"),
305 _title_phi_2("NA"),
306 _title_phi_3("NA"),
307 _title_phi_12("NA"),
308 _title_phi_13("NA"),
309 _title_phi_23("NA"),
310 _title_phi_123("NA"),
311
312 vsZ("NA"),
313 vsM("NA"),
314 vsPt("NA"),
315 vsPhi("NA"), 
316 vsEta("NA"), 
317 vsEtaPhi("NA")
318 {
319   printf("Default constructor called \n");
320   
321   printf("passed \n ");
322   
323 }
324
325 AliAnalysisTask3PCorrelations::AliAnalysisTask3PCorrelations(const TString & name)
326 : AliAnalysisTaskSE(name),
327 fAODEvent(0), 
328 fESDEvent(0),             //! ESD Event 
329 fInputHandler(0),
330 _outputHistoList(0),
331 _twoPi         ( 2.0 * 3.1415927),
332 _eventCount    ( 0), 
333 _debugLevel    ( 0),
334 _singlesOnly   ( 0), 
335 _useWeights    ( 0), 
336 _rejectPileup  ( 1), 
337 _rejectPairConversion ( 0), 
338 _vertexZMin           ( -10), 
339 _vertexZMax           (  10), 
340 _vertexXYMin          ( -10),
341 _vertexXYMax          (  10),
342 _centralityMethod     (  4),
343 _centralityMin        (  0.),
344 _centralityMax        (  0.),
345 _dcaZMin              ( -3), 
346 _dcaZMax              (  3.), 
347 _dcaXYMin             ( -3.), 
348 _dcaXYMax             (  3.),
349 _dedxMin              ( 0),
350 _dedxMax              ( 100000),
351 _nClusterMin          ( 70), 
352 _trackFilterBit       ( 128),
353 _field    ( 1.),
354 _nTracks  ( 0 ),
355 _mult0    ( 0 ),
356 _mult1    ( 0 ),
357 _mult2    ( 0 ),
358 _mult3    ( 0 ),
359 _mult4    ( 0 ),
360 _mult5    ( 0 ),
361 _mult6    ( 0 ),
362 arraySize ( 2000),
363 _id_1(0),       
364 _charge_1(0),    
365 _iPhi_1(0),    
366 _pt_1(0),       
367 _px_1(0),      
368 _py_1(0),      
369 _pz_1(0),      
370 _correction_1(0),   
371 _dedx_1(0),   
372 _id_2(0),      
373 _charge_2(0),    
374 _iPhi_2(0),    
375 _pt_2(0),      
376 _px_2(0),      
377 _py_2(0),      
378 _pz_2(0),      
379 _correction_2(0),   
380 _dedx_2(0),   
381 _id_3(0),      
382 _charge_3(0),    
383 _iPhi_3(0),    
384 _pt_3(0),      
385 _px_3(0),      
386 _py_3(0),      
387 _pz_3(0),      
388 _correction_3(0),   
389 _dedx_3(0),   
390
391 _correctionWeight_1P(0),   
392 _correctionWeight_1M(0),   
393 _correctionWeight_2P(0),   
394 _correctionWeight_2M(0),   
395 _correctionWeight_3P(0),   
396 _correctionWeight_3M(0),   
397 _nBins_M0(500),       _min_M0(0),        _max_M0(10000),          _width_M0(20),
398 _nBins_M1(500),       _min_M1(0),        _max_M1(10000),          _width_M1(20),
399 _nBins_M2(500),       _min_M2(0),        _max_M2(10000),          _width_M2(20),
400 _nBins_M3(500),       _min_M3(0),        _max_M3(10000),          _width_M3(20),
401 _nBins_M4(100),       _min_M4(0),        _max_M4(1),              _width_M4(0.01),
402 _nBins_M5(100),       _min_M5(0),        _max_M5(1),              _width_M5(0.01),
403 _nBins_M6(100),       _min_M6(0),        _max_M6(1),              _width_M6(0.01),
404 _nBins_vsM(0),        _min_vsM(0),       _max_vsM(0),             _width_vsM(0),
405
406 _nBins_vertexZ(40),   _min_vertexZ(-10), _max_vertexZ(10),        _width_vertexZ(0.5),
407
408 _nBins_pt_1(8),       _min_pt_1(2.),    _max_pt_1(10.0),          _width_pt_1(1.0),
409 _nBins_phi_1(72),     _min_phi_1(0),     _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
410 _nBins_eta_1(1),      _min_eta_1(-0.25), _max_eta_1(0.25),        _width_eta_1(0.5),
411 _nBins_etaPhi_1(0), 
412 _nBins_etaPhiPt_1(0),
413 _nBins_zEtaPhiPt_1(0),
414
415 _nBins_pt_2(5),       _min_pt_2(0.5),     _max_pt_2(1.0),          _width_pt_2(0.1),
416 _nBins_phi_2(72),     _min_phi_2(0),      _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
417 _nBins_eta_2(1),      _min_eta_2(-1),     _max_eta_2(-0.7),        _width_eta_2(0.3),
418 _nBins_etaPhi_2(0), 
419 _nBins_etaPhiPt_2(0),
420 _nBins_zEtaPhiPt_2(0),
421
422 _nBins_pt_3(5),       _min_pt_3(0.5),     _max_pt_3(1.0),          _width_pt_3(0.1),
423 _nBins_phi_3(72),     _min_phi_3(0),      _max_phi_3(2.*3.1415927),_width_phi_3(2.*3.1415927/72),
424 _nBins_eta_3(1),      _min_eta_3(0.7),    _max_eta_3(1.),          _width_eta_3(0.3),
425 _nBins_etaPhi_3(0), 
426 _nBins_etaPhiPt_3(0),
427 _nBins_zEtaPhiPt_3(0),
428
429 _nBins_phi_12(0),
430 _nBins_phi_13(0),
431 _nBins_phi_23(0),
432 _nBins_phi_123(0),
433
434 _nBins_etaPhi_12(0),
435 _nBins_etaPhi_13(0),
436 _nBins_etaPhi_23(0),
437 _nBins_etaPhi_123(0),
438
439 __n1_1(0),
440 __n1_2(0),
441 __n1_3(0),
442 __n2_12(0),   
443 __n2_13(0),   
444 __n2_23(0),   
445 __n3_123(0),   
446 __n1Nw_1(0),
447 __n1Nw_2(0),
448 __n1Nw_3(0),
449 __n2Nw_12(0),   
450 __n2Nw_13(0),   
451 __n2Nw_23(0),   
452 __n3Nw_123(0),  
453
454 __n1_1_vsPt(0),
455 __n1_1_vsPhi(0), 
456 __n1_1P_vsZEtaPhiPt(0),
457 __n1_1M_vsZEtaPhiPt(0),
458
459 __n1_2_vsPt(0),
460 __n1_2_vsPhi(0), 
461 __n1_2P_vsZEtaPhiPt(0),
462 __n1_2M_vsZEtaPhiPt(0),
463
464 __n1_3_vsPt(0),
465 __n1_3_vsPhi(0), 
466 __n1_3P_vsZEtaPhiPt(0),
467 __n1_3M_vsZEtaPhiPt(0),
468
469 __n2_12_vsPhi(0),
470 __n2_13_vsPhi(0),
471 __n2_23_vsPhi(0),
472 __n3_123_vsPhi(0),
473
474 _weight_1P      ( 0    ),
475 _weight_1M      ( 0    ),
476 _weight_2P      ( 0    ),
477 _weight_2M      ( 0    ),
478 _weight_3P      ( 0    ),
479 _weight_3M      ( 0    ),
480 _eventAccounting ( 0),
481 _m0 ( 0),
482 _m1 ( 0),
483 _m2 ( 0),
484 _m3 ( 0),
485 _m4 ( 0),
486 _m5 ( 0),
487 _m6 ( 0),
488 _vertexZ        ( 0),
489
490 _n1_1_vsPhi     ( 0),
491 _n1_1_vsM          ( 0), // w/ weight
492 _n1Nw_1_vsM        ( 0), // w/o weight
493
494 _n1_1_vsPt      ( 0),         
495 _n1_1P_vsZVsEtaVsPhiVsPt ( 0),
496 _n1_1M_vsZVsEtaVsPhiVsPt ( 0),
497 _dedxVsP_1         ( 0),
498 _corrDedxVsP_1     ( 0),
499 _betaVsP_1         ( 0),
500
501 _n1_2_vsPhi        ( 0),
502 _n1_2_vsM          ( 0),
503 _n1Nw_2_vsM        ( 0),
504 _n1_2_vsPt         ( 0),       
505 _n1_2P_vsZVsEtaVsPhiVsPt ( 0), 
506 _n1_2M_vsZVsEtaVsPhiVsPt ( 0), 
507 _dedxVsP_2         ( 0),
508 _corrDedxVsP_2     ( 0),
509 _betaVsP_2         ( 0),
510
511 _n1_3_vsPhi        ( 0),
512 _n1_3_vsM          ( 0),
513 _n1Nw_3_vsM        ( 0),
514 _n1_3_vsPt         ( 0),       
515 _n1_3P_vsZVsEtaVsPhiVsPt ( 0), 
516 _n1_3M_vsZVsEtaVsPhiVsPt ( 0), 
517 _dedxVsP_3         ( 0),
518 _corrDedxVsP_3     ( 0),
519 _betaVsP_3         ( 0),
520
521 _n2_12_vsPhi    ( 0),  
522 _n2_13_vsPhi    ( 0),  
523 _n2_23_vsPhi    ( 0),  
524 _n3_123_vsPhi   ( 0),  
525
526 _n2_12_vsM         ( 0),        
527 _n2_13_vsM         ( 0),        
528 _n2_23_vsM         ( 0),        
529 _n3_123_vsM        ( 0),        
530
531 _n2Nw_12_vsM       ( 0),        
532 _n2Nw_13_vsM       ( 0),        
533 _n2Nw_23_vsM       ( 0),        
534 _n3Nw_123_vsM      ( 0),        
535
536 n1Name("NA"),
537 n2Name("NA"),
538 n3Name("NA"),
539 n1NwName("NA"),
540 n2NwName("NA"),
541 n3NwName("NA"),
542 n1n1Name("NA"),
543 n1n1n1Name("NA"),
544 n2n1Name("NA"),
545 r1Name("NA"),
546 r2Name("NA"),
547 r3Name("NA"),
548 r2r1Name("NA"),
549 c2Name("NA"),
550 c3Name("NA"),
551 d3Name("NA"),
552 p3Name("NA"),
553
554 intR2Name("NA"),
555 binCorrName("NA"),
556 intBinCorrName("NA"),
557
558 countsName("NA"),
559 part_1_Name("NA"),
560 part_1P_Name("NA"),
561 part_1M_Name("NA"),
562 part_2_Name("NA"),
563 part_2P_Name("NA"),
564 part_2M_Name("NA"),
565 part_3_Name("NA"),
566 part_3P_Name("NA"),
567 part_3M_Name("NA"),
568 pair_12_Name("NA"),
569 pair_13_Name("NA"),
570 pair_23_Name("NA"),
571 triplet_Name("NA"),
572
573 _title_counts("NA"),
574 _title_m0("NA"),
575 _title_m1("NA"),
576 _title_m2("NA"),
577 _title_m3("NA"),
578 _title_m4("NA"),
579 _title_m5("NA"),
580 _title_m6("NA"),
581
582 _title_phi_1("NA"),
583 _title_phi_2("NA"),
584 _title_phi_3("NA"),
585
586 _title_phi_12("NA"),
587 _title_phi_13("NA"),
588 _title_phi_23("NA"),
589 _title_phi_123("NA"),
590
591 vsZ("NA"),
592 vsM("NA"),
593 vsPt("NA"),
594 vsPhi("NA"), 
595 vsEta("NA"), 
596 vsEtaPhi("NA")
597 {
598   printf("2nd constructor called ");
599   
600   DefineOutput(0, TList::Class());
601   
602   printf("passed  ");
603   
604 }
605
606 AliAnalysisTask3PCorrelations::~AliAnalysisTask3PCorrelations()
607 {
608   // no ops --- at least for now....
609 }
610
611 void AliAnalysisTask3PCorrelations::UserCreateOutputObjects()
612 {
613   cout<< "AliAnalysisTask3PCorrelations::CreateOutputObjects() Starting " << endl;
614   OpenFile(0);
615   _outputHistoList = new TList();
616   _outputHistoList->SetOwner();
617   
618   //if (_useWeights) DefineInput(2, TList::Class());
619   
620   //Setup the parameters of the histograms
621   _nBins_M0 = 500; _min_M0   = 0.;    _max_M0    = 5000.;  _width_M0 = (_max_M0-_min_M0)/_nBins_M0;
622   _nBins_M1 = 500; _min_M1   = 0.;    _max_M1    = 5000.;  _width_M1 = (_max_M1-_min_M1)/_nBins_M1;
623   _nBins_M2 = 500; _min_M2   = 0.;    _max_M2    = 5000.;  _width_M2 = (_max_M2-_min_M2)/_nBins_M2;
624   _nBins_M3 = 500; _min_M3   = 0.;    _max_M3    = 5000.;  _width_M3 = (_max_M3-_min_M3)/_nBins_M3;
625   _nBins_M4 = 100; _min_M4   = 0.;    _max_M4    = 100.;   _width_M4 = (_max_M4-_min_M4)/_nBins_M4;
626   _nBins_M5 = 100; _min_M5   = 0.;    _max_M5    = 100.;   _width_M5 = (_max_M5-_min_M5)/_nBins_M5;
627   _nBins_M6 = 100; _min_M6   = 0.;    _max_M6    = 100.;   _width_M6 = (_max_M6-_min_M6)/_nBins_M6;
628   
629   _min_vertexZ       = _vertexZMin;  
630   _max_vertexZ       = _vertexZMax;  
631   _width_vertexZ     = 0.5;
632   _nBins_vertexZ     = int(0.5+ (_max_vertexZ - _min_vertexZ)/_width_vertexZ); 
633   
634   _nBins_pt_1        = int(0.5+ (_max_pt_1 -_min_pt_1 )/_width_pt_1); 
635   _nBins_eta_1       = int(0.5+ (_max_eta_1-_min_eta_1)/_width_eta_1);  
636   _width_phi_1       = (_max_phi_1  - _min_phi_1)  /_nBins_phi_1;
637   _nBins_etaPhi_1    = _nBins_phi_1    * _nBins_eta_1;
638   _nBins_zEtaPhiPt_1 = _nBins_vertexZ  * _nBins_phi_1 * _nBins_eta_1 * _nBins_pt_1;
639   
640   _nBins_pt_2        =  int(0.5+ (_max_pt_2 -_min_pt_2 )/_width_pt_2);  
641   _nBins_eta_2       = int(0.5+ (_max_eta_2-_min_eta_2)/_width_eta_2); 
642   _width_phi_2       = (_max_phi_2  - _min_phi_2)  /_nBins_phi_2;
643   _nBins_etaPhi_2    = _nBins_phi_2    * _nBins_eta_2;
644   _nBins_zEtaPhiPt_2 = _nBins_vertexZ   * _nBins_phi_2 * _nBins_eta_2 * _nBins_pt_2;
645   
646   _nBins_pt_3        =  int(0.5+ (_max_pt_3 -_min_pt_3 )/_width_pt_3);  
647   _nBins_eta_3       = int(0.5+ (_max_eta_3-_min_eta_3)/_width_eta_3); 
648   _width_phi_3       = (_max_phi_3  - _min_phi_3)  /_nBins_phi_3;
649   _nBins_etaPhi_3    = _nBins_phi_3    * _nBins_eta_3;
650   _nBins_zEtaPhiPt_3 = _nBins_vertexZ   * _nBins_phi_3 * _nBins_eta_3 * _nBins_pt_3;
651   
652   
653   _nBins_etaPhi_12   = _nBins_etaPhi_1 * _nBins_etaPhi_2;
654   _nBins_etaPhi_13   = _nBins_etaPhi_1 * _nBins_etaPhi_3;
655   _nBins_etaPhi_23   = _nBins_etaPhi_2 * _nBins_etaPhi_3;
656   _nBins_etaPhi_123  = _nBins_etaPhi_1 * _nBins_etaPhi_2 * _nBins_etaPhi_3;
657   
658   _nBins_phi_12   = _nBins_phi_1 * _nBins_phi_2;
659   _nBins_phi_13   = _nBins_phi_1 * _nBins_phi_3;
660   _nBins_phi_23   = _nBins_phi_2 * _nBins_phi_3;
661   _nBins_phi_123  = _nBins_phi_1 * _nBins_phi_2 * _nBins_phi_3;
662   
663   
664   //setup the work arrays
665   if (_singlesOnly)
666     {
667     __n1_1_vsPt         = getDoubleArray(_nBins_pt_1,         0.);
668     __n1_1P_vsZEtaPhiPt = getFloatArray (_nBins_zEtaPhiPt_1,  0.);
669     __n1_1M_vsZEtaPhiPt = getFloatArray (_nBins_zEtaPhiPt_1,  0.);
670     __n1_2_vsPt         = getDoubleArray(_nBins_pt_2,         0.);
671     __n1_2P_vsZEtaPhiPt = getFloatArray (_nBins_zEtaPhiPt_2,  0.);
672     __n1_2M_vsZEtaPhiPt = getFloatArray (_nBins_zEtaPhiPt_2,  0.);
673     __n1_3_vsPt         = getDoubleArray(_nBins_pt_3,         0.);
674     __n1_3P_vsZEtaPhiPt = getFloatArray (_nBins_zEtaPhiPt_3,  0.);
675     __n1_3M_vsZEtaPhiPt = getFloatArray (_nBins_zEtaPhiPt_3,  0.);
676     }
677   else
678     {
679     _id_1          = new int[arraySize];   
680     //_charge_1      = new int[arraySize]; 
681     _iPhi_1        = new int[arraySize]; 
682     //_pt_1          = new float[arraySize];    
683     //_px_1          = new float[arraySize];   
684     //_py_1          = new float[arraySize];   
685     //_pz_1          = new float[arraySize];   
686     _correction_1  = new float[arraySize];
687     //_dedx_1        = new float[arraySize];
688     __n1_1_vsPhi   = getDoubleArray(_nBins_phi_1,    0.);
689     
690     _id_2          = new int[arraySize];   
691     //_charge_2      = new int[arraySize]; 
692     _iPhi_2        = new int[arraySize]; 
693     //_pt_2          = new float[arraySize];   
694     //_px_2          = new float[arraySize];   
695     //_py_2          = new float[arraySize];   
696     //_pz_2          = new float[arraySize];   
697     _correction_2  = new float[arraySize];
698     //_dedx_2        = new float[arraySize];
699     __n1_2_vsPhi   = getDoubleArray(_nBins_phi_2,       0.);
700     
701     _id_3          = new int[arraySize];   
702     //_charge_3      = new int[arraySize]; 
703     _iPhi_3        = new int[arraySize]; 
704     //_pt_3          = new float[arraySize];   
705     //_px_3          = new float[arraySize];   
706     //_py_3          = new float[arraySize];   
707     //_pz_3          = new float[arraySize];   
708     _correction_3  = new float[arraySize];
709     //_dedx_3        = new float[arraySize];
710     __n1_3_vsPhi        = getDoubleArray(_nBins_phi_3,       0.);
711     
712     __n2_12_vsPhi    = getFloatArray(_nBins_phi_12,    0.);
713     __n2_13_vsPhi    = getFloatArray(_nBins_phi_13,    0.);
714     __n2_23_vsPhi    = getFloatArray(_nBins_phi_23,    0.);
715     __n3_123_vsPhi   = getFloatArray(_nBins_phi_123,   0.);
716     }
717   
718   
719   // Setup all the labels needed.
720   part_1_Name   = "_1";
721   part_1P_Name  = "_1P";
722   part_1M_Name  = "_1M";
723   part_2_Name   = "_2";
724   part_2P_Name  = "_2P";
725   part_2M_Name  = "_2M";
726   part_3_Name   = "_3";
727   part_3P_Name  = "_3P";
728   part_3M_Name  = "_3M";
729   pair_12_Name  = "_12";
730   pair_13_Name  = "_13";
731   pair_23_Name  = "_23";
732   triplet_Name   = "_123";
733   
734   n1Name     = "n1";
735   n2Name     = "n2";
736   n3Name     = "n3";
737   n1NwName   = "n1Nw";
738   n2NwName   = "n2Nw";
739   n3NwName   = "n3Nw";
740   
741   r1Name     = "r1";
742   r2Name     = "r2";
743   r3Name     = "r3";
744   r2r1Name   = "r2r1";
745   c2Name     = "c2";
746   c3Name     = "c3";
747   d3Name     = "d3";
748   p3Name     = "p3";
749   
750   intR2Name       = "intR2";
751   binCorrName     = "binCorr";
752   intBinCorrName  = "intBinCorr";
753   
754   
755   _title_counts = "yield";
756   _title_m0     = "M_{0}";
757   _title_m1     = "M_{1}";
758   _title_m2     = "M_{2}";
759   _title_m3     = "M_{3}";
760   _title_m4     = "V0Centrality";
761   _title_m5     = "TrkCentrality";
762   _title_m6     = "SpdCentrality";
763   
764   _title_phi_1       = "#varphi_{1} (radian)";
765   _title_phi_2       = "#varphi_{2} (radian)";
766   _title_phi_3       = "#varphi_{3} (radian)";
767   _title_phi_12   = "#eta_{1}#times#varphi_{1}#times#eta_{2}#times#varphi_{2}";
768   _title_phi_13   = "#eta_{1}#times#varphi_{1}#times#eta_{3}#times#varphi_{3}";
769   _title_phi_23   = "#eta_{2}#times#varphi_{2}#times#eta_{3}#times#varphi_{3}";
770   
771   vsZ         = "_vsZ";
772   vsM         = "_vsM";
773   vsPt         = "_vsPt";
774   vsPhi        = "_vsPhi"; 
775   vsEta        = "_vsEta"; 
776   vsEtaPhi     = "vsPhi"; 
777   
778   
779   if (_useWeights)
780     {
781     int iZ, iPhi, iPt;
782     int iZ1,iPhi1,iPt1;
783     int a, b;
784     if (_weight_1P)
785       {
786       _correctionWeight_1P = new float[_nBins_vertexZ*_nBins_etaPhi_1*_nBins_pt_1];
787       _correctionWeight_1M = new float[_nBins_vertexZ*_nBins_etaPhi_1*_nBins_pt_1];
788       a = _nBins_pt_1;
789       b = _nBins_etaPhi_1*_nBins_pt_1;
790       for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
791         {
792         for (iPhi=0,iPhi1=1; iPhi<_nBins_etaPhi_1; iPhi++, iPhi1++)
793           {
794           for (iPt=0,iPt1=1; iPt<_nBins_pt_1; iPt++, iPt1++)
795             {
796             _correctionWeight_1P[iZ*b+iPhi*a+iPt] = _weight_1P->GetBinContent(iZ1,iPhi1,iPt1);
797             _correctionWeight_1M[iZ*b+iPhi*a+iPt] = _weight_1M->GetBinContent(iZ1,iPhi1,iPt1);
798             }      
799           }
800         }
801       } // _weight_1
802     else
803       {
804       AliError("AliAnalysisTask3PCorrelations:: _weight_1 is a null pointer.");
805       return;
806       }
807     if (_weight_2P)
808       {
809       _correctionWeight_2P = new float[_nBins_vertexZ*_nBins_etaPhi_2*_nBins_pt_2];
810       _correctionWeight_2M = new float[_nBins_vertexZ*_nBins_etaPhi_2*_nBins_pt_2];
811       a = _nBins_pt_2;
812       b = _nBins_etaPhi_2*_nBins_pt_2;
813       for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
814         {
815         for (iPhi=0,iPhi1=1; iPhi<_nBins_etaPhi_2; iPhi++, iPhi1++)
816           {
817           for (iPt=0,iPt1=1; iPt<_nBins_pt_2; iPt++, iPt1++)
818             {
819             _correctionWeight_2P[iZ*b+iPhi*a+iPt] = _weight_2P->GetBinContent(iZ1,iPhi1,iPt1);
820             _correctionWeight_2M[iZ*b+iPhi*a+iPt] = _weight_2M->GetBinContent(iZ1,iPhi1,iPt1);
821             }      
822           }
823         }
824       } // _weight_2
825     else
826       {
827       AliError("AliAnalysisTask3PCorrelations:: _weight_2 is a null pointer.");
828       return;
829       }
830     if (_weight_3P)
831       {
832       _correctionWeight_3P = new float[_nBins_vertexZ*_nBins_etaPhi_3*_nBins_pt_3];
833       _correctionWeight_3M = new float[_nBins_vertexZ*_nBins_etaPhi_3*_nBins_pt_3];
834       a = _nBins_pt_3;
835       b = _nBins_etaPhi_3*_nBins_pt_3;
836       for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
837         {
838         for (iPhi=0,iPhi1=1; iPhi<_nBins_etaPhi_2; iPhi++, iPhi1++)
839           {
840           for (iPt=0,iPt1=1; iPt<_nBins_pt_3; iPt++, iPt1++)
841             {
842             _correctionWeight_3P[iZ*b+iPhi*a+iPt] = _weight_3P->GetBinContent(iZ1,iPhi1,iPt1);
843             _correctionWeight_3M[iZ*b+iPhi*a+iPt] = _weight_3M->GetBinContent(iZ1,iPhi1,iPt1);
844             }      
845           }
846         }
847       } // _weight_2
848     else
849       {
850       AliError("AliAnalysisTask3PCorrelations:: _weight_2 is a null pointer.");
851       return;
852       }
853     }
854   
855   createHistograms();
856   PostData(0,_outputHistoList);
857   
858   cout<< "AliAnalysisTask3PCorrelations::CreateOutputObjects() DONE " << endl;
859   
860 }
861
862 void  AliAnalysisTask3PCorrelations::createHistograms()
863 {
864   AliInfo(" AliAnalysisTask3PCorrelations::createHistoHistograms() Creating Event Histos");
865   TString name;
866   
867   // bin index : what it is...
868   //         0 :  number of event submitted
869   //         1 :  number accepted by centrality cut
870   //         2 :  number accepted by centrality cut and z cut
871   //         3 :  total number of particles that satisfied filter 1
872   //         4 :  total number of particles that satisfied filter 2
873   //         5 :  total number of particles that satisfied filter 3
874   name = "eventAccounting";
875   _eventAccounting      = createHisto1D(name,name,10, -0.5, 9.5, "event Code", _title_counts);
876   
877   name = "m0"; _m0      = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m0, _title_counts);
878   name = "m1"; _m1      = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m1, _title_counts);
879   name = "m2"; _m2      = createHisto1D(name,name,_nBins_M2, _min_M2, _max_M2, _title_m2, _title_counts);
880   name = "m3"; _m3      = createHisto1D(name,name,_nBins_M3, _min_M3, _max_M3, _title_m3, _title_counts);
881   name = "m4"; _m4      = createHisto1D(name,name,_nBins_M4, _min_M4, _max_M4, _title_m4, _title_counts);
882   name = "m5"; _m5      = createHisto1D(name,name,_nBins_M5, _min_M5, _max_M5, _title_m5, _title_counts);
883   name = "m6"; _m6      = createHisto1D(name,name,_nBins_M6, _min_M6, _max_M6, _title_m6, _title_counts);
884   name = "zV"; _vertexZ = createHisto1D(name,name,_nBins_vertexZ, _min_vertexZ, _max_vertexZ, "z-Vertex (cm)", _title_counts);
885   
886   TString title_vsM;
887   //_centralityMethod
888   switch (_centralityMethod)
889   {
890     case 0: _nBins_vsM = _nBins_M0; _min_vsM = _min_M0; _max_vsM = _max_M0; _width_vsM = _width_M0; title_vsM = _title_m0; break;
891     case 1: _nBins_vsM = _nBins_M1; _min_vsM = _min_M1; _max_vsM = _max_M1; _width_vsM = _width_M1; title_vsM = _title_m1; break;
892     case 2: _nBins_vsM = _nBins_M2; _min_vsM = _min_M2; _max_vsM = _max_M2; _width_vsM = _width_M2; title_vsM = _title_m2; break;
893     case 3: _nBins_vsM = _nBins_M3; _min_vsM = _min_M3; _max_vsM = _max_M3; _width_vsM = _width_M3; title_vsM = _title_m3; break;
894     case 4: _nBins_vsM = _nBins_M4; _min_vsM = _min_M4; _max_vsM = _max_M4; _width_vsM = _width_M4; title_vsM = _title_m4; break;
895     case 5: _nBins_vsM = _nBins_M5; _min_vsM = _min_M5; _max_vsM = _max_M5; _width_vsM = _width_M5; title_vsM = _title_m5; break;
896     case 6: _nBins_vsM = _nBins_M6; _min_vsM = _min_M6; _max_vsM = _max_M6; _width_vsM = _width_M6; title_vsM = _title_m6; break;
897   }
898   
899   AliInfo(" AliAnalysisTask3PCorrelations::createHistoHistograms() Creating Part 1 Histos");
900   
901   if (_singlesOnly) 
902     {
903     name = n1Name+part_1_Name+vsPt;               _n1_1_vsPt               = createHisto1F(name,name, _nBins_pt_1,  _min_pt_1,  _max_pt_1,   "pt1", "counts");
904     name = "dedxVsP_1";                           _dedxVsP_1               = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
905     name = "betaVsP_1";                           _betaVsP_1               = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
906     name = n1Name+part_1P_Name+vsZ+vsEtaPhi+vsPt; _n1_1P_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_1, 0., double(_nBins_etaPhi_1), _nBins_pt_1, _min_pt_1, _max_pt_1, "zVertex", _title_phi_1,  "pt1");
907     name = n1Name+part_1M_Name+vsZ+vsEtaPhi+vsPt; _n1_1M_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_1, 0., double(_nBins_etaPhi_1), _nBins_pt_1, _min_pt_1, _max_pt_1, "zVertex", _title_phi_1,  "pt1");
908     
909     name = n1Name+part_2_Name+vsPt;               _n1_2_vsPt               = createHisto1F(name,name, _nBins_pt_2,  _min_pt_2,  _max_pt_2,   "pt2",  "counts");
910     name = "dedxVsP_2";                           _dedxVsP_2               = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
911     name = "betaVsP_2";                           _betaVsP_2               = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
912     name = n1Name+part_2P_Name+vsZ+vsEtaPhi+vsPt; _n1_2P_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_2, 0., double(_nBins_etaPhi_2), _nBins_pt_2, _min_pt_2, _max_pt_2, "zVertex", _title_phi_2,  "pt2");
913     name = n1Name+part_2M_Name+vsZ+vsEtaPhi+vsPt; _n1_2M_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_2, 0., double(_nBins_etaPhi_2), _nBins_pt_2, _min_pt_2, _max_pt_2, "zVertex", _title_phi_2,  "pt2");
914     
915     name = n1Name+part_3_Name+vsPt;               _n1_3_vsPt               = createHisto1F(name,name, _nBins_pt_3,  _min_pt_3,  _max_pt_3,   "pt3",  "counts");
916     name = "dedxVsP_3";                           _dedxVsP_3               = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
917     name = "betaVsP_3";                           _betaVsP_3               = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
918     name = n1Name+part_3P_Name+vsZ+vsEtaPhi+vsPt; _n1_3P_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_3, 0., double(_nBins_etaPhi_3), _nBins_pt_3, _min_pt_3, _max_pt_3, "zVertex", _title_phi_3,  "pt3");
919     name = n1Name+part_3M_Name+vsZ+vsEtaPhi+vsPt; _n1_3M_vsZVsEtaVsPhiVsPt = createHisto3F(name,name, _nBins_vertexZ,_min_vertexZ,_max_vertexZ, _nBins_etaPhi_3, 0., double(_nBins_etaPhi_3), _nBins_pt_3, _min_pt_3, _max_pt_3, "zVertex", _title_phi_3,  "pt3");
920     }
921   else
922     {
923     name = n1Name+part_1_Name+vsPhi;  _n1_1_vsPhi   = createHisto1F(name,name, _nBins_phi_1, _min_phi_1, _max_phi_1, _title_phi_1,  "counts");
924     name = n1Name+part_2_Name+vsPhi;  _n1_2_vsPhi   = createHisto1F(name,name, _nBins_phi_2, _min_phi_2, _max_phi_2, _title_phi_2,  "counts");
925     name = n1Name+part_3_Name+vsPhi;  _n1_3_vsPhi   = createHisto1F(name,name, _nBins_phi_3, _min_phi_3, _max_phi_3, _title_phi_3,  "counts");
926     
927     name = n1Name+part_1_Name + vsM;  _n1_1_vsM     = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{1}>");
928     name = n1Name+part_2_Name + vsM;  _n1_2_vsM     = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{2}>");
929     name = n1Name+part_3_Name + vsM;  _n1_3_vsM     = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{3}>");
930     
931     name = n1NwName+part_1_Name+vsM;  _n1Nw_1_vsM   = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{1}>");
932     name = n1NwName+part_2_Name+vsM;  _n1Nw_2_vsM   = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{2}>");
933     name = n1NwName+part_3_Name+vsM;  _n1Nw_3_vsM   = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{3}>");
934     
935     name = n2Name+pair_12_Name+vsPhi; _n2_12_vsPhi  = createHisto2F(name,name, _nBins_phi_1, _min_phi_1, _max_phi_1, _nBins_phi_2, _min_phi_2, _max_phi_2, _title_phi_1, _title_phi_2, "<N_{12}>");
936     name = n2Name+pair_13_Name+vsPhi; _n2_13_vsPhi  = createHisto2F(name,name, _nBins_phi_1, _min_phi_1, _max_phi_1, _nBins_phi_3, _min_phi_3, _max_phi_3, _title_phi_2, _title_phi_3, "<N_{13}>");        
937     name = n2Name+pair_23_Name+vsPhi; _n2_23_vsPhi  = createHisto2F(name,name, _nBins_phi_2, _min_phi_2, _max_phi_2, _nBins_phi_3, _min_phi_3, _max_phi_3, _title_phi_1, _title_phi_3, "<N_{23}>");
938     name = n3Name+triplet_Name+vsPhi; _n3_123_vsPhi = createHisto3F(name,name, _nBins_phi_1, _min_phi_1, _max_phi_1, _nBins_phi_2, _min_phi_2, _max_phi_2, _nBins_phi_3, _min_phi_3, _max_phi_3, _title_phi_1, _title_phi_2, _title_phi_3);        
939     
940     name = n2Name+pair_12_Name+vsM;  _n2_12_vsM     = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{12}>");
941     name = n2Name+pair_13_Name+vsM;  _n2_13_vsM     = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{13}>");
942     name = n2Name+pair_23_Name+vsM;  _n2_23_vsM     = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{23}>");
943     name = n3Name+pair_23_Name+vsM;  _n3_123_vsM    = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{123}>");
944     name = n2NwName+pair_12_Name+vsM; _n2Nw_12_vsM  = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{12}>");
945     name = n2NwName+pair_13_Name+vsM; _n2Nw_13_vsM  = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{13}>");
946     name = n2NwName+pair_23_Name+vsM; _n2Nw_23_vsM  = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{23}>");
947     name = n3NwName+triplet_Name+vsM; _n3Nw_123_vsM = createProfile(name,name, _nBins_vsM,   _min_vsM,   _max_vsM,  title_vsM,      "<N_{123}>");
948     
949     }
950   
951   AliInfo(" AliAnalysisTask3PCorrelations::createHistoHistograms() All Done"); 
952 }
953 //-----------------------//
954
955 void  AliAnalysisTask3PCorrelations::finalizeHistograms()
956 {
957   
958   AliInfo("AliAnalysisTask3PCorrelations::finalizeHistograms() starting");
959   AliInfo(Form("CorrelationAnalyzers::finalizeHistograms()   _eventCount : %d",int(_eventCount)));
960   
961   if (_singlesOnly) 
962     {
963     fillHistoWithArray(_n1_1_vsPt,               __n1_1_vsPt,         _nBins_pt_1);
964     fillHistoWithArray(_n1_1P_vsZVsEtaVsPhiVsPt, __n1_1P_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
965     fillHistoWithArray(_n1_1M_vsZVsEtaVsPhiVsPt, __n1_1M_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
966     
967     fillHistoWithArray(_n1_2_vsPt,               __n1_2_vsPt,         _nBins_pt_2);
968     fillHistoWithArray(_n1_2P_vsZVsEtaVsPhiVsPt, __n1_2P_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_2, _nBins_pt_2);
969     fillHistoWithArray(_n1_2M_vsZVsEtaVsPhiVsPt, __n1_2M_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_2, _nBins_pt_2);
970     
971     fillHistoWithArray(_n1_3_vsPt,               __n1_3_vsPt,         _nBins_pt_3);
972     fillHistoWithArray(_n1_3P_vsZVsEtaVsPhiVsPt, __n1_3P_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_3, _nBins_pt_3);
973     fillHistoWithArray(_n1_3M_vsZVsEtaVsPhiVsPt, __n1_3M_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_3, _nBins_pt_3);
974     
975     }
976   else
977     {
978     fillHistoWithArray(_n1_1_vsPhi,    __n1_1_vsPhi,   _nBins_phi_1);
979     fillHistoWithArray(_n1_2_vsPhi,    __n1_2_vsPhi,   _nBins_phi_2);
980     fillHistoWithArray(_n1_3_vsPhi,    __n1_3_vsPhi,   _nBins_phi_3);
981     fillHistoWithArray(_n2_12_vsPhi,   __n2_12_vsPhi,  _nBins_phi_1, _nBins_phi_2);
982     fillHistoWithArray(_n2_13_vsPhi,   __n2_13_vsPhi,  _nBins_phi_1, _nBins_phi_3);
983     fillHistoWithArray(_n2_23_vsPhi,   __n2_23_vsPhi,  _nBins_phi_2, _nBins_phi_3);
984     fillHistoWithArray(_n3_123_vsPhi,  __n3_123_vsPhi, _nBins_phi_1, _nBins_phi_2, _nBins_phi_3);
985     }
986   AliInfo("AliAnalysisTask3PCorrelations::finalizeHistograms()  Done ");
987 }
988 //--------------//
989
990
991 void  AliAnalysisTask3PCorrelations::UserExec(Option_t */*option*/)
992 {
993   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - Starting!!!!" << endl;
994   
995   int    k1,k2,k3;
996   int    iPhi, iEta, iEtaPhi, iPt, charge;
997   float  q, p, phi, pt, eta, corr, dedx, px, py, pz;
998   int    ij_12, ij_13, ij_23, ijk_123;
999   int    id_1, id_2, id_3, iPhi_1, iPhi_2, iPhi_3;
1000   float  corr_1, corr_2, corr_3;
1001   float  corr_12, corr_13, corr_23, corr_123;
1002   int    iVertex, iVertexP1, iVertexP2, iVertexP3;
1003   int    iZEtaPhiPt;
1004   //float  massElecSq = 2.5e-7;
1005   double b[2];
1006   double bCov[3];
1007   const  AliAODVertex*  vertex;
1008   int    nClus;
1009   bool   bitOK;
1010   
1011   // count all events looked at here
1012   _eventCount++;
1013   
1014   if (_eventAccounting)
1015     {
1016     _eventAccounting->Fill(0);// count all calls to this function
1017     }
1018   else
1019     {
1020     cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - !_eventAccounting" << endl;
1021     cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - Skip this task" << endl;
1022     return;
1023     }
1024   
1025   //Get pointer to current event
1026   //fESDEvent = dynamic_cast<AliESDEvent*> (InputEvent());
1027   fAODEvent = dynamic_cast<AliAODEvent*> (InputEvent());
1028   
1029   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 2" << endl;
1030   
1031   if(!fAODEvent)
1032     {
1033     AliError("AliAnalysisTask3PCorrelations::Exec(Option_t *option) !fAODEvent"); 
1034     return;
1035     }
1036   
1037   
1038   float v0Centr  = -999.;
1039   float trkCentr = -999.;
1040   float spdCentr = -999.;
1041   float vertexX  = -999;
1042   float vertexY  = -999;
1043   float vertexZ  = -999;
1044   float vertexXY = -999;
1045   float dcaZ     = -999;
1046   float dcaXY    = -999;
1047   float centrality = -999;
1048   k1 = k2 = k3 = 0;
1049   __n1_1 = __n1_2 = __n1_3 = __n1Nw_1 = __n1Nw_2 = __n1Nw_3 = 0;
1050   
1051   _eventAccounting->Fill(1);// count all calls to this function with a valid pointer
1052   
1053   //Centrality
1054   AliCentrality* centralityObject =  fAODEvent->GetHeader()->GetCentralityP();
1055   if (centralityObject)
1056     {
1057     // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 6" << endl;
1058     
1059     v0Centr  = centralityObject->GetCentralityPercentile("V0M");
1060     trkCentr = centralityObject->GetCentralityPercentile("TRK"); 
1061     spdCentr = centralityObject->GetCentralityPercentile("CL1");
1062     // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 7" << endl;
1063     
1064     }
1065   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 8" << endl;
1066   
1067   _nTracks  = fAODEvent->GetNTracks();
1068   _mult0    = 0;
1069   _mult1    = 0;
1070   _mult2    = 0;
1071   _mult3    = _nTracks; 
1072   _mult4    = v0Centr;
1073   _mult5    = trkCentr;
1074   _mult6    = spdCentr;
1075   _field    = fAODEvent->GetMagneticField(); 
1076   
1077   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - _field:" << _field << endl;
1078   
1079   //_centralityMethod
1080   switch (_centralityMethod)
1081   {
1082     case 0: centrality = _mult0; break;
1083     case 1: centrality = _mult1; break;
1084     case 2: centrality = _mult2; break;
1085     case 3: centrality = _mult3; break;
1086     case 4: centrality = _mult4; break;
1087     case 5: centrality = _mult5; break;
1088     case 6: centrality = _mult6; break;
1089   }
1090   
1091   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 10" << endl;
1092   
1093   //filter on centrality
1094   if ( centrality < _centralityMin ||
1095       centrality > _centralityMax)
1096     {
1097     // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 11" << endl;
1098     
1099     // we dont analyze this event here... 
1100     return;
1101     }
1102   _eventAccounting->Fill(2);// count all events with right centrality
1103                             // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 12" << endl;
1104   
1105   // filter on z and xy vertex
1106   vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertexSPD();
1107   if (!vertex || vertex->GetNContributors()<1)
1108     {
1109     vertexZ  = -999;
1110     vertexXY = -999;
1111     // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - No valid vertex object or poor vertex" << endl;
1112     }
1113   else
1114     { 
1115       vertexX = vertex->GetX();
1116       vertexY = vertex->GetY();
1117       vertexZ = vertex->GetZ();
1118       vertexXY = sqrt(vertexX*vertexX+vertexY*vertexY);
1119     }
1120   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 13" << endl;
1121   
1122   // cout << "vertexZ : " << vertexZ  << endl;
1123   // cout << "vertexXY: " << vertexXY << endl;
1124   // cout << "_vertexZMin: " << _vertexZMin << endl;
1125   // cout << "_vertexZMax: " << _vertexZMax << endl;
1126   // cout << "_vertexXYMin: " << _vertexXYMin << endl;
1127   // cout << "_vertexXYMax: " << _vertexXYMax << endl;
1128   
1129   if (!vertex ||
1130       vertexZ  < _vertexZMin  || 
1131       vertexZ  > _vertexZMax  ||
1132       vertexXY < _vertexXYMin || 
1133       vertexXY > _vertexXYMax)
1134     return;
1135   
1136   
1137   // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 14" << endl;
1138   
1139   iVertex = int((vertexZ-_min_vertexZ)/_width_vertexZ);
1140   iVertexP1 = iVertex*_nBins_etaPhiPt_1;
1141   iVertexP2 = iVertex*_nBins_etaPhiPt_2;
1142   iVertexP3 = iVertex*_nBins_etaPhiPt_3;
1143   if (iVertex<0 || iVertex>=_nBins_vertexZ)
1144     {
1145     AliError("AliAnalysisTask3PCorrelations::Exec(Option_t *option) iVertex<0 || iVertex>=_nBins_vertexZ ");
1146     return;
1147     }
1148   _eventAccounting->Fill(3);// count all calls to this function with a valid pointer
1149                             // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 15" << endl;
1150   
1151   
1152   // loop over all particles for singles
1153   //debug() << "_nTracks:"<< _nTracks << endl;
1154   //reset single particle counters
1155   for (int iTrack=0; iTrack< _nTracks; iTrack++)
1156     {
1157     // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - 16 iTrack:" << iTrack << endl;
1158     
1159     AliAODTrack * t = (AliAODTrack *) fAODEvent->GetTrack(iTrack);
1160     if (!t) 
1161       {
1162       AliError(Form("AliAnalysisTask3PCorrelations::Exec(Option_t *option) No track ofr iTrack=%d", iTrack));
1163       continue;
1164       }
1165     
1166     q      = t->Charge();
1167     charge = int(q);
1168     phi    = t->Phi();
1169     p      = t->P();
1170     pt     = t->Pt();
1171     px     = t->Px();
1172     py     = t->Py();
1173     pz     = t->Pz();
1174     eta    = t->Eta();
1175     dedx   = t->GetTPCsignal();
1176     nClus  = t->GetTPCNcls();
1177     bitOK  = t->TestFilterBit(_trackFilterBit);
1178     if (!bitOK ||
1179         dedx  < _dedxMin   ||
1180         dedx  > _dedxMax   ||
1181         nClus < _nClusterMin) continue;
1182     
1183     // cout << "_trackFilterBit:" << _trackFilterBit << " Track returns:" << bitOK << endl;
1184     // cout << "   pt:" << pt   << " _min_pt_1:" << _min_pt_1 << " _max_pt_1:" << _max_pt_1<< endl;
1185     // cout << "  phi:" << phi  << endl;
1186     // cout << "  eta:" << eta  << " _min_eta_1:" << _min_eta_1 << " _max_eta_1:" << _max_eta_1<< endl;
1187     // cout << " dedx:" << dedx << " _dedxMin:" << _dedxMin << " _dedxMax:" << _dedxMax << endl;
1188     // cout << "nclus:" << nClus<< " _nClusterMin:" << _nClusterMin << endl;
1189     if (pt       >=  _min_pt_1  && 
1190         pt       <   _max_pt_1  && 
1191         eta      >=  _min_eta_1 && 
1192         eta      <   _max_eta_1) 
1193       {
1194       // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - check vertex  for 1:" << endl;
1195       // Get the dca information
1196       if (t->PropagateToDCA(vertex, _field, 100., b, bCov) )
1197         {
1198         dcaXY = b[0];
1199         dcaZ  = b[1];
1200         } 
1201       else
1202         {
1203         dcaXY = -999999;
1204         dcaZ  = -999999;
1205         }
1206       // cout << "1 dcaZ:" << dcaZ << " _dcaZMin:" << _dcaZMin << " _dcaZMax:" << _dcaZMax << endl;
1207       // cout << "1 dcaXY:" << dcaXY << " _dcaXYMin:" << _dcaXYMin << " _dcaXYMax:" << _dcaXYMax << endl;
1208       
1209       // skip track if DCA too large
1210       if (dcaZ     >=  _dcaZMin && 
1211           dcaZ     <   _dcaZMax && 
1212           dcaXY    >=  _dcaXYMin && 
1213           dcaXY    <   _dcaXYMax)
1214         continue; //track does not have a valid DCA
1215                   // cout << "keep track:" << endl;
1216       
1217       iPhi   = int( phi/_width_phi_1);
1218       // cout << " AliAnalysisTask3PCorrelations::analyze(Event * event) -1- iTrack:" << iTrack<< endl<< "pt:" << pt << " phi:" <<  phi << " eta:" << eta << endl;
1219       if (iPhi<0 || iPhi>=_nBins_phi_1 ) 
1220         {
1221         AliWarning("AliAnalysisTask3PCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1222         return;
1223         }
1224       
1225       iEta    = int((eta-_min_eta_1)/_width_eta_1); 
1226       if (iEta<0 || iEta>=_nBins_eta_1) 
1227         {
1228         AliWarning(Form("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1229         continue;
1230         }
1231       iPt     = int((pt -_min_pt_1 )/_width_pt_1 ); 
1232       if (iPt<0  || iPt >=_nBins_pt_1)
1233         {
1234         AliWarning(Form("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1235         continue;
1236         }
1237       iEtaPhi    = iEta*_nBins_phi_1+iPhi;
1238       iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
1239       
1240       if (charge>0 && _correctionWeight_1P)
1241         corr = _correctionWeight_1P[iZEtaPhiPt];
1242       else if (charge<0 && _correctionWeight_1M)
1243         corr = _correctionWeight_1M[iZEtaPhiPt];
1244       else
1245         corr = 1;
1246       // cout << "all good; process track:" << endl;
1247       if (_singlesOnly)
1248         {
1249         __n1_1_vsPt[iPt]                += corr;
1250         if (charge>0) 
1251           __n1_1P_vsZEtaPhiPt[iZEtaPhiPt] += corr;
1252         else if (charge>0) 
1253           __n1_1M_vsZEtaPhiPt[iZEtaPhiPt] += corr;
1254         }
1255       else
1256         {
1257         _id_1[k1]            = iTrack;     
1258         //_charge_1[k1]        = charge;
1259         _iPhi_1[k1]          = iPhi; 
1260         _correction_1[k1]    = corr; 
1261         __n1_1               += corr;
1262         __n1_1_vsPhi[iPhi]   += corr; 
1263         __n1Nw_1             += 1;
1264         ++k1;
1265         if (k1>=arraySize)
1266           {
1267           AliError(Form("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
1268           return;
1269           }
1270         }
1271       }
1272     
1273     else if (pt  >= _min_pt_2  && 
1274              pt  <  _max_pt_2  && 
1275              eta >= _min_eta_2 && 
1276              eta <  _max_eta_2)  
1277       {
1278       // Get the dca information
1279       // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - check vertex for 2:" << endl;
1280       if (t->PropagateToDCA(vertex, _field, 100., b, bCov) )
1281         {
1282         dcaXY = b[0];
1283         dcaZ  = b[1];
1284         } 
1285       else
1286         {
1287         dcaXY = -999999;
1288         dcaZ  = -999999;
1289         }
1290       // cout << "2 dcaZ:" << dcaZ << " _dcaZMin:" << _dcaZMin << " _dcaZMax:" << _dcaZMax << endl;
1291       // cout << "2 dcaXY:" << dcaXY << " _dcaXYMin:" << _dcaXYMin << " _dcaXYMax:" << _dcaXYMax << endl;
1292       // skip track if DCA too large
1293       if (dcaZ     >=  _dcaZMin && 
1294           dcaZ     <   _dcaZMax && 
1295           dcaXY    >=  _dcaXYMin && 
1296           dcaXY    <   _dcaXYMax)
1297         continue; //track does not have a valid DCA
1298       
1299       iPhi   = int( phi/_width_phi_2);
1300       //cout << " AliAnalysisTask3PCorrelations::analyze(Event * event) -1- iTrack:" << iTrack  << endl
1301       //<< "pt:" << pt << " phi:" <<  phi << " eta:" << eta << endl;
1302       if (iPhi<0 || iPhi>=_nBins_phi_2 ) 
1303         {
1304         AliWarning("AliAnalysisTask3PCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1305         return;
1306         }
1307       
1308       iEta    = int((eta-_min_eta_2)/_width_eta_2);
1309       if (iEta<0 || iEta>=_nBins_eta_2) 
1310         {
1311         AliWarning(Form("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1312         continue;
1313         }
1314       iPt     = int((pt -_min_pt_2 )/_width_pt_2 ); 
1315       if (iPt<0  || iPt >=_nBins_pt_2)
1316         {
1317         AliWarning(Form("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1318         continue;
1319         }
1320       
1321       iEtaPhi    = iEta*_nBins_phi_2+iPhi;
1322       iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
1323       if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
1324         {
1325         AliWarning("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
1326         continue;
1327         }
1328       
1329       if (charge>0 && _correctionWeight_2P)
1330         corr = _correctionWeight_2P[iZEtaPhiPt];
1331       else if (charge<0 && _correctionWeight_2M)
1332         corr = _correctionWeight_2M[iZEtaPhiPt];
1333       else
1334         corr = 1;
1335       
1336       if (_singlesOnly)
1337         {
1338         //_dedxVsP_2->Fill(p*q,dedx);
1339         __n1_2_vsPt[iPt]               += corr;          // cout << "step 15" << endl;
1340         if (charge>0) 
1341           __n1_2P_vsZEtaPhiPt[iZEtaPhiPt] += corr;
1342         else if (charge>0) 
1343           __n1_2M_vsZEtaPhiPt[iZEtaPhiPt] += corr;          }
1344       else
1345         {
1346         _id_2[k2]           = iTrack;         // cout << "step 1" << endl;
1347                                               //_charge_2[k2]       = charge;         // cout << "step 2" << endl;
1348         _iPhi_2[k2]         = iPhi;        // cout << "step 3" << endl;
1349         _correction_2[k2]   = corr;           // cout << "step 9" << endl;
1350         __n1_2              += corr;          // cout << "step 10" << endl;
1351         __n1Nw_2            += 1;
1352         __n1_2_vsPhi[iPhi]  += corr;          // cout << "step 11" << endl;
1353         ++k2;
1354         if (k2>=arraySize)
1355           {
1356           AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize)); 
1357           return;
1358           }
1359         }
1360       }
1361     
1362     else if (pt  >= _min_pt_3  && 
1363              pt  <  _max_pt_3  && 
1364              eta >= _min_eta_3 && 
1365              eta <  _max_eta_3)  
1366       {
1367       // Get the dca information
1368       // cout << "AliAnalysisTask3PCorrelations::UserExec(Option_t *option) - check vertex for 2:" << endl;
1369       if (t->PropagateToDCA(vertex, _field, 100., b, bCov) )
1370         {
1371         dcaXY = b[0];
1372         dcaZ  = b[1];
1373         } 
1374       else
1375         {
1376         dcaXY = -999999;
1377         dcaZ  = -999999;
1378         }
1379       // skip track if DCA too large
1380       if (dcaZ     >=  _dcaZMin && 
1381           dcaZ     <   _dcaZMax && 
1382           dcaXY    >=  _dcaXYMin && 
1383           dcaXY    <   _dcaXYMax)
1384         continue; //track does not have a valid DCA
1385       
1386       iPhi   = int( phi/_width_phi_3);
1387       if (iPhi<0 || iPhi>=_nBins_phi_3 ) 
1388         {
1389         AliWarning("AliAnalysisTask3PCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_3");
1390         return;
1391         }
1392       
1393       iEta    = int((eta-_min_eta_3)/_width_eta_3);
1394       if (iEta<0 || iEta>=_nBins_eta_3) 
1395         {
1396         AliWarning(Form("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1397         continue;
1398         }
1399       iEtaPhi    = iEta*_nBins_phi_3+iPhi;
1400       iZEtaPhiPt = iVertexP3 + iEtaPhi*_nBins_pt_3 + iPt;
1401       if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_3)
1402         {
1403         AliWarning("AliAnalysisTask3PCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_3");
1404         continue;
1405         }
1406       if (charge>0 && _correctionWeight_3P)
1407         corr = _correctionWeight_3P[iZEtaPhiPt];
1408       else if (charge<0 && _correctionWeight_3M)
1409         corr = _correctionWeight_3M[iZEtaPhiPt];
1410       else
1411         corr = 1;
1412       
1413       if (_singlesOnly)
1414         {
1415         __n1_3_vsPt[iPt]               += corr;          // cout << "step 15" << endl;
1416         if (charge>0) 
1417           __n1_3P_vsZEtaPhiPt[iZEtaPhiPt] += corr;
1418         else if (charge>0) 
1419           __n1_3M_vsZEtaPhiPt[iZEtaPhiPt] += corr;          }
1420       else
1421         {
1422         _id_3[k3]            = iTrack;         // cout << "step 1" << endl;
1423                                                //_charge_3[k3]        = charge;         // cout << "step 2" << endl;
1424         _iPhi_3[k3]          = iPhi;        // cout << "step 3" << endl;
1425         _correction_3[k3]    = corr;           // cout << "step 9" << endl;
1426         __n1_3               += corr;          // cout << "step 10" << endl;
1427         __n1Nw_3             += 1;
1428         __n1_3_vsPhi[iPhi]   += corr;          // cout << "step 11" << endl;
1429         ++k3;
1430         if (k3>=arraySize)
1431           {
1432           AliWarning(Form("-W- k3 >=arraySize; arraySize: %d",arraySize)); 
1433           return;
1434           }
1435         }
1436       }
1437     
1438     // cout << "done with track" << endl;
1439     } //iTrack
1440   
1441   // cout << "Filling histograms now" << endl;
1442   _m0->Fill(_mult0);
1443   _m1->Fill(_mult1);
1444   _m2->Fill(_mult2);
1445   _m3->Fill(_mult3);
1446   _m4->Fill(_mult4);
1447   _m5->Fill(_mult5);
1448   _m6->Fill(_mult6);
1449   _vertexZ->Fill(vertexZ);
1450   
1451   
1452   //if singels only selected, do not fill pair histograms.
1453   if (_singlesOnly) 
1454     {
1455     
1456     }
1457   else
1458     {
1459     // reset pair counters
1460     _n1_1_vsM->Fill(    centrality, __n1_1);
1461     _n1Nw_1_vsM->Fill(  centrality, __n1Nw_1);
1462     _n1_2_vsM->Fill(    centrality, __n1_2);
1463     _n1Nw_2_vsM->Fill(  centrality, __n1Nw_2);
1464     _n1_3_vsM->Fill(    centrality, __n1_3);
1465     _n1Nw_3_vsM->Fill(  centrality, __n1Nw_3);
1466
1467     __n2_12 = __n2Nw_12 = __n2_13 = __n2Nw_13 = __n2_23 = __n2Nw_23 = __n3_123 = __n3Nw_123 = 0;
1468     
1469     for (int i1=0; i1<k1; i1++)
1470       {
1471       // cout << "         i1:" << i1 << endl;
1472       id_1    = _id_1[i1]; 
1473       iPhi_1  = _iPhi_1[i1];
1474       corr_1  = _correction_1[i1];
1475       //1 and 2
1476       for (int i2=0;i2<k2;i2++)
1477         {        
1478           id_2 = _id_2[i2];
1479           if (id_1!=id_2)
1480             {
1481             iPhi_2  = _iPhi_2[i2];  
1482             corr_2  = _correction_2[i2];
1483             corr_12 = corr_1*corr_2;
1484             ij_12 = iPhi_1*_nBins_phi_1 + iPhi_2;
1485             __n2_12_vsPhi[ij_12] += corr_12;
1486             __n2_12              += corr_12;
1487             __n2Nw_12            += 1;
1488             for (int i3=0;i3<k3;i3++)
1489               {        
1490                 id_3 = _id_3[i3];
1491                 if (id_1!=id_3 && id_2!=id_3)
1492                   {
1493                   iPhi_3   = _iPhi_3[i3];  
1494                   corr_3   = _correction_3[i3];
1495                   corr_123 = corr_12*corr_3;
1496                   ijk_123  = ij_12 *_nBins_phi_2 + iPhi_3;
1497                   __n3_123_vsPhi[ijk_123] += corr_123;
1498                   __n3_123                += corr_123;
1499                   __n3Nw_123              += 1;
1500                   }
1501               } //i3       
1502             }
1503         } //i2       
1504           //1 and 3
1505       for (int i3=0;i3<k3;i3++)
1506         {        
1507           id_3 = _id_3[i3];
1508           if (id_1!=id_3)
1509             {
1510             iPhi_3  = _iPhi_3[i3];  
1511             corr_3  = _correction_3[i3];
1512             corr_13 = corr_1*corr_3;
1513             ij_13   = iPhi_1*_nBins_phi_1 + iPhi_3;
1514             __n2_13_vsPhi[ij_13] += corr_13;
1515             __n2_13              += corr_13;
1516             __n2Nw_13            += 1;
1517             }
1518         } //i3       
1519       } //i1  
1520     
1521     //2 & 3
1522     for (int i2=0; i2<k2; i2++)
1523       {
1524       id_2    = _id_2[i2]; 
1525       iPhi_2  = _iPhi_2[i2];
1526       corr_2  = _correction_2[i2];
1527       for (int i3=0;i3<k3;i3++)
1528         {        
1529           id_3 = _id_3[i3];
1530           if (id_2!=id_3)
1531             {
1532             iPhi_3  = _iPhi_3[i3];  
1533             corr_3  = _correction_3[i3];
1534             corr_23 = corr_2*corr_3;
1535             ij_23   = iPhi_2*_nBins_phi_2 + iPhi_3;
1536             __n2_23_vsPhi[ij_23] += corr_23;
1537             __n2_23              += corr_23;
1538             __n2Nw_23            += 1;
1539             }
1540         } //i2
1541       }
1542     _n2_12_vsM->Fill(   centrality, __n2_12);
1543     _n2_13_vsM->Fill(   centrality, __n2_13);
1544     _n2_23_vsM->Fill(   centrality, __n2_23);
1545     _n3_123_vsM->Fill(  centrality, __n3_123);
1546     _n2Nw_12_vsM->Fill( centrality, __n2Nw_12);
1547     _n2Nw_13_vsM->Fill( centrality, __n2Nw_13);
1548     _n2Nw_23_vsM->Fill( centrality, __n2Nw_23);
1549     _n3Nw_123_vsM->Fill(centrality, __n3_123);
1550     }
1551   
1552   // cout << "Event Done " << endl;
1553   PostData(0,_outputHistoList);
1554   
1555 }
1556
1557 void   AliAnalysisTask3PCorrelations::FinishTaskOutput()
1558 {
1559   cout << "AliAnalysisTask3PCorrelations::FinishTaskOutput() Starting." << endl;
1560   finalizeHistograms();
1561   PostData(0,_outputHistoList);
1562   cout << "AliAnalysisTask3PCorrelations::FinishTaskOutput() Done." << endl;
1563 }
1564
1565 void   AliAnalysisTask3PCorrelations::Terminate(Option_t* /*option*/)
1566 {
1567   // no ops
1568 }
1569
1570
1571 //Tools
1572 //===================================================================================================
1573 void  AliAnalysisTask3PCorrelations::fillHistoWithArray(TH1 * h, float * array, int size)
1574 {
1575   int i, i1;
1576   float v1, ev1, v2, ev2, sum, esum;
1577   for (i=0, i1=1; i<size; ++i,++i1)
1578     {
1579     v1  = array[i]; ev1 = sqrt(v1);
1580     v2  = h->GetBinContent(i1);
1581     ev2 = h->GetBinError(i1);
1582     sum = v1 + v2;
1583     esum = sqrt(ev1*ev1+ev2*ev2);
1584     h->SetBinContent(i1,sum);
1585     h->SetBinError(i1,esum);
1586     }
1587 }
1588
1589 void  AliAnalysisTask3PCorrelations::fillHistoWithArray(TH2 * h, float * array, int size1, int size2)
1590 {
1591   int i, i1;
1592   int j, j1;
1593   float v1, ev1, v2, ev2, sum, esum;
1594   for (i=0, i1=1; i<size1; ++i,++i1)
1595     {
1596     for (j=0, j1=1; j<size2; ++j,++j1)
1597       {
1598       v1  = array[i*size2+j]; ev1 = sqrt(v1);
1599       v2  = h->GetBinContent(i1,j1);
1600       ev2 = h->GetBinError(i1,j1);
1601       sum = v1 + v2;
1602       esum = sqrt(ev1*ev1+ev2*ev2);
1603       h->SetBinContent(i1,j1,sum);
1604       h->SetBinError(i1,j1,esum);
1605       }
1606     }
1607 }
1608
1609 void  AliAnalysisTask3PCorrelations::fillHistoWithArray(TH3 * h, float * array, int size1, int size2, int size3)
1610 {
1611   int i, i1;
1612   int j, j1;
1613   int k, k1;
1614   float v1, ev1, v2, ev2, sum, esum;
1615   int size23 = size2*size3;
1616   for (i=0, i1=1; i<size1; ++i,++i1)
1617     {
1618     for (j=0, j1=1; j<size2; ++j,++j1)
1619       {
1620       for (k=0, k1=1; k<size3; ++k,++k1)
1621         {
1622         v1  = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
1623         v2  = h->GetBinContent(i1,j1,k1);
1624         ev2 = h->GetBinError(i1,j1,k1);
1625         sum = v1 + v2;
1626         esum = sqrt(ev1*ev1+ev2*ev2);
1627         h->SetBinContent(i1,j1,k1,sum);
1628         h->SetBinError(i1,j1,k1,esum);
1629         }
1630       }
1631     }
1632 }
1633
1634 void  AliAnalysisTask3PCorrelations::fillHistoWithArray(TH1 * h, double * array, int size)
1635 {
1636   int i, i1;
1637   double v1, ev1, v2, ev2, sum, esum;
1638   for (i=0, i1=1; i<size; ++i,++i1)
1639     {
1640     v1  = array[i]; ev1 = sqrt(v1);
1641     v2  = h->GetBinContent(i1);
1642     ev2 = h->GetBinError(i1);
1643     sum = v1 + v2;
1644     esum = sqrt(ev1*ev1+ev2*ev2);
1645     h->SetBinContent(i1,sum);
1646     h->SetBinError(i1,esum);
1647     }
1648 }
1649
1650 void  AliAnalysisTask3PCorrelations::fillHistoWithArray(TH2 * h, double * array, int size1, int size2)
1651 {
1652   int i, i1;
1653   int j, j1;
1654   double v1, ev1, v2, ev2, sum, esum;
1655   for (i=0, i1=1; i<size1; ++i,++i1)
1656     {
1657     for (j=0, j1=1; j<size2; ++j,++j1)
1658       {
1659       v1  = array[i*size2+j]; ev1 = sqrt(v1);
1660       v2  = h->GetBinContent(i1,j1);
1661       ev2 = h->GetBinError(i1,j1);
1662       sum = v1 + v2;
1663       esum = sqrt(ev1*ev1+ev2*ev2);
1664       h->SetBinContent(i1,j1,sum);
1665       h->SetBinError(i1,j1,esum);
1666       }
1667     }
1668 }
1669
1670 void  AliAnalysisTask3PCorrelations::fillHistoWithArray(TH3 * h, double * array, int size1, int size2, int size3)
1671 {
1672   int i, i1;
1673   int j, j1;
1674   int k, k1;
1675   double v1, ev1, v2, ev2, sum, esum;
1676   int size23 = size2*size3;
1677   for (i=0, i1=1; i<size1; ++i,++i1)
1678     {
1679     for (j=0, j1=1; j<size2; ++j,++j1)
1680       {
1681       for (k=0, k1=1; k<size3; ++k,++k1)
1682         {
1683         v1  = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
1684         v2  = h->GetBinContent(i1,j1,k1);
1685         ev2 = h->GetBinError(i1,j1,k1);
1686         sum = v1 + v2;
1687         esum = sqrt(ev1*ev1+ev2*ev2);
1688         h->SetBinContent(i1,j1,k1,sum);
1689         h->SetBinError(i1,j1,k1,esum);
1690         }
1691       }
1692     }
1693 }
1694
1695 //________________________________________________________________________
1696 double *  AliAnalysisTask3PCorrelations::getDoubleArray(int size, double v)
1697 {
1698   /// Allocate an array of type double with n values
1699   /// Initialize the array to the given value
1700   double * array = new double [size];
1701   for (int i=0;i<size;++i) array[i]=v;
1702   return array;
1703 }
1704
1705 //________________________________________________________________________
1706 float *  AliAnalysisTask3PCorrelations::getFloatArray(int size, float v)
1707 {
1708   /// Allocate an array of type float with n values
1709   /// Initialize the array to the given value
1710   float * array = new float [size];
1711   for (int i=0;i<size;++i) array[i]=v;
1712   return array;
1713 }
1714
1715
1716 //________________________________________________________________________
1717 TH1D * AliAnalysisTask3PCorrelations::createHisto1D(const TString &  name, const TString &  title, 
1718                                                     int n, double xMin, double xMax, 
1719                                                     const TString &  xTitle, const TString &  yTitle)
1720 {
1721   //CreateHisto new 1D historgram
1722   AliInfo(Form("createHisto 1D histo %s   nBins: %d  xMin: %f   xMax: %f",name.Data(),n,xMin,xMax));
1723   TH1D * h = new TH1D(name,title,n,xMin,xMax);
1724   h->GetXaxis()->SetTitle(xTitle);
1725   h->GetYaxis()->SetTitle(yTitle);
1726   addToList(h);
1727   return h;
1728 }
1729
1730
1731 //________________________________________________________________________
1732 TH1D * AliAnalysisTask3PCorrelations::createHisto1D(const TString &  name, const TString &  title, 
1733                                                     int n, double * bins, 
1734                                                     const TString &  xTitle, const TString &  yTitle)
1735 {
1736   AliInfo(Form("createHisto 1D histo %s   with %d non uniform nBins",name.Data(),n));
1737   TH1D * h = new TH1D(name,title,n,bins);
1738   h->GetXaxis()->SetTitle(xTitle);
1739   h->GetYaxis()->SetTitle(yTitle);
1740   addToList(h);
1741   return h;
1742 }
1743
1744
1745 //________________________________________________________________________
1746 TH2D * AliAnalysisTask3PCorrelations::createHisto2D(const TString &  name, const TString &  title, 
1747                                                     int nx, double xMin, double xMax, int ny, double yMin, double yMax, 
1748                                                     const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
1749 {
1750   AliInfo(Form("createHisto 2D histo %s  nx: %d  xMin: %f10.4 xMax: %f10.4  ny: %d   yMin: %f10.4 yMax: %f10.4",name.Data(),nx,xMin,xMax,ny,yMin,yMax));
1751   TH2D * h = new TH2D(name,title,nx,xMin,xMax,ny,yMin,yMax);
1752   h->GetXaxis()->SetTitle(xTitle);
1753   h->GetYaxis()->SetTitle(yTitle);
1754   h->GetZaxis()->SetTitle(zTitle);
1755   addToList(h);
1756   return h;
1757 }
1758
1759 //________________________________________________________________________
1760 TH2D * AliAnalysisTask3PCorrelations::createHisto2D(const TString &  name, const TString &  title, 
1761                                                     int nx, double* xbins, int ny, double yMin, double yMax, 
1762                                                     const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
1763 {
1764   AliInfo(Form("createHisto 2D histo %s   with %d non uniform nBins",name.Data(),nx));
1765   TH2D * h;
1766   h = new TH2D(name,title,nx,xbins,ny,yMin,yMax);
1767   h->GetXaxis()->SetTitle(xTitle);
1768   h->GetYaxis()->SetTitle(yTitle);
1769   h->GetZaxis()->SetTitle(zTitle);
1770   addToList(h);
1771   return h;
1772 }
1773
1774 //// F /////
1775 //________________________________________________________________________
1776 TH1F * AliAnalysisTask3PCorrelations::createHisto1F(const TString &  name, const TString &  title, 
1777                                                     int n, double xMin, double xMax, 
1778                                                     const TString &  xTitle, const TString &  yTitle)
1779 {
1780   //CreateHisto new 1D historgram
1781   AliInfo(Form("createHisto 1D histo %s   nBins: %d  xMin: %f   xMax: %f",name.Data(),n,xMin,xMax));
1782   TH1F * h = new TH1F(name,title,n,xMin,xMax);
1783   h->GetXaxis()->SetTitle(xTitle);
1784   h->GetYaxis()->SetTitle(yTitle);
1785   addToList(h);
1786   return h;
1787 }
1788
1789
1790 //________________________________________________________________________
1791 TH1F * AliAnalysisTask3PCorrelations::createHisto1F(const TString &  name, const TString &  title, 
1792                                                     int n, double * bins, 
1793                                                     const TString &  xTitle, const TString &  yTitle)
1794 {
1795   AliInfo(Form("createHisto 1D histo %s   with %d non uniform nBins",name.Data(),n));
1796   TH1F * h = new TH1F(name,title,n,bins);
1797   h->GetXaxis()->SetTitle(xTitle);
1798   h->GetYaxis()->SetTitle(yTitle);
1799   addToList(h);
1800   return h;
1801 }
1802
1803
1804 //________________________________________________________________________
1805 TH2F * AliAnalysisTask3PCorrelations::createHisto2F(const TString &  name, const TString &  title, 
1806                                                     int nx, double xMin, double xMax, int ny, double yMin, double yMax, 
1807                                                     const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
1808 {
1809   AliInfo(Form("createHisto 2D histo %s  nx: %d  xMin: %f10.4 xMax: %f10.4  ny: %d   yMin: %f10.4 yMax: %f10.4",name.Data(),nx,xMin,xMax,ny,yMin,yMax));
1810   TH2F * h = new TH2F(name,title,nx,xMin,xMax,ny,yMin,yMax);
1811   h->GetXaxis()->SetTitle(xTitle);
1812   h->GetYaxis()->SetTitle(yTitle);
1813   h->GetZaxis()->SetTitle(zTitle);
1814   addToList(h);
1815   return h;
1816 }
1817
1818 //________________________________________________________________________
1819 TH2F * AliAnalysisTask3PCorrelations::createHisto2F(const TString &  name, const TString &  title, 
1820                                                     int nx, double* xbins, int ny, double yMin, double yMax, 
1821                                                     const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
1822 {
1823   AliInfo(Form("createHisto 2D histo %s   with %d non uniform nBins",name.Data(),nx));
1824   TH2F * h;
1825   h = new TH2F(name,title,nx,xbins,ny,yMin,yMax);
1826   h->GetXaxis()->SetTitle(xTitle);
1827   h->GetYaxis()->SetTitle(yTitle);
1828   h->GetZaxis()->SetTitle(zTitle);
1829   addToList(h);
1830   return h;
1831 }
1832
1833
1834 //________________________________________________________________________
1835 TH3F * AliAnalysisTask3PCorrelations::createHisto3F(const TString &  name, const TString &  title, 
1836                                                     int nx, double xMin, double xMax, 
1837                                                     int ny, double yMin, double yMax, 
1838                                                     int nz, double zMin, double zMax, 
1839                                                     const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
1840 {
1841   AliInfo(Form("createHisto 2D histo %s  nx: %d  xMin: %f10.4 xMax: %f10.4  ny: %d   yMin: %f10.4 yMax: %f10.4 nz: %d   zMin: %f10.4 zMax: %f10.4",name.Data(),nx,xMin,xMax,ny,yMin,yMax,nz,zMin,zMax));
1842   TH3F * h = new TH3F(name,title,nx,xMin,xMax,ny,yMin,yMax,nz,zMin,zMax);
1843   h->GetXaxis()->SetTitle(xTitle);
1844   h->GetYaxis()->SetTitle(yTitle);
1845   h->GetZaxis()->SetTitle(zTitle);
1846   addToList(h);
1847   return h;
1848 }
1849
1850 //________________________________________________________________________
1851 TProfile * AliAnalysisTask3PCorrelations::createProfile(const TString & name, const TString & description, 
1852                                                         int nx,double xMin,double xMax,
1853                                                         const TString &  xTitle, const TString &  yTitle)
1854 {
1855   AliInfo(Form("createHisto 1D profile %s   nBins: %d  xMin: %f10.4 xMax: %f10.4",name.Data(),nx,xMin,xMax));
1856   TProfile * h = new TProfile(name,description,nx,xMin,xMax);
1857   h->GetXaxis()->SetTitle(xTitle);
1858   h->GetYaxis()->SetTitle(yTitle);
1859   addToList(h);
1860   return h;  
1861 }
1862
1863 //________________________________________________________________________
1864 TProfile * AliAnalysisTask3PCorrelations::createProfile(const TString &  name,const TString &  description, 
1865                                                         int nx,  double* bins,
1866                                                         const TString &  xTitle, const TString &  yTitle)
1867 {
1868   AliInfo(Form("createHisto 1D profile %s  with %d non uniform bins",name.Data(),nx));
1869   TProfile * h = new TProfile(name,description,nx,bins);
1870   h->GetXaxis()->SetTitle(xTitle);
1871   h->GetYaxis()->SetTitle(yTitle);
1872   addToList(h);
1873   return h;
1874 }
1875
1876
1877 void   AliAnalysisTask3PCorrelations::addToList(TH1 *h)
1878 {
1879   if (_outputHistoList)
1880     {
1881     _outputHistoList->Add(h);
1882     }
1883   else
1884     cout << "addToList(TH1 *h) _outputHistoList is null!!!!! Shoudl abort ship" << endl;
1885   
1886 }
1887
1888
1889