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