PWGCF/Correlations/DPhi/AliAnalysisTaskDptDptCorrelations.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 "TChain.h"
  19 #include "TList.h"
  20 #include "TFile.h"
  21 #include "TTree.h"
  22 #include "TH1D.h"
  23 #include "TH2D.h"
  24 #include "TH3D.h"
  25 #include "THnSparse.h"
  26 #include "TCanvas.h"
  27
  28 #include <TROOT.h>
  29 #include <TChain.h>
  30 #include <TFile.h>
  31 #include <TList.h>
  32 #include <TMath.h>
  33 #include <TTree.h>
  34 #include <TH1F.h>
  35 #include <TH2F.h>
  36 #include <TH3F.h>
  37 #include <TProfile.h>
  38 #include <TH1D.h>
  39 #include <TH2D.h>
  40 #include <TH3D.h>
  41 #include <TRandom.h>
  42 #include "AliAnalysisManager.h"
  43
  44 #include "AliAODHandler.h"
  45 #include "AliAODInputHandler.h"
  46 #include "AliInputEventHandler.h"
  47 #include "AliLog.h"
  48 #include "AliESDEvent.h"
  49 #include "AliESDInputHandler.h"
  50 #include "AliMultiplicity.h"
  51 #include "AliCentrality.h"
  52 #include "AliAnalysisTaskDptDptCorrelations.h"
  53
  54 #include "AliPID.h"
  55 #include "AliPIDResponse.h"
  56
  57 #include "AliESDVertex.h"
  58 #include "AliESDEvent.h"
  59 #include "AliESDInputHandler.h"
  60 #include "AliAODEvent.h"
  61 #include "AliAODTrack.h"
  62 #include "AliAODInputHandler.h"
  63 #include "AliESD.h"
  64 #include "AliESDEvent.h"
  65 #include "AliAODEvent.h"
  66 #include "AliStack.h"
  67 #include "AliESDtrackCuts.h"
  68 #include "AliAODMCHeader.h"
  69
  70
  71 #include "AliGenHijingEventHeader.h"
  72 #include "AliGenEventHeader.h"
  73 #include "AliPID.h"
  74 #include "AliAODPid.h"
  75 #include "AliPIDResponse.h"
  76 #include "AliAODpidUtil.h"
  77 #include "AliPIDCombined.h"
  78
  79
  80 ClassImp(AliAnalysisTaskDptDptCorrelations)
  81
  82 AliAnalysisTaskDptDptCorrelations::AliAnalysisTaskDptDptCorrelations()
  83 : AliAnalysisTaskSE(),
  84 fAODEvent(0),
  85 fESDEvent(0),             //! ESD Event
  86 fInputHandler(0),
  87 fPIDResponse(0x0),
  88 _outputHistoList(0),
  89 _twoPi         ( 2.0 * 3.1415927),
  90 _eventCount    ( 0),
  91 _debugLevel    ( 0),
  92 _singlesOnly   ( 0),
  93 _useWeights    ( 0),
  94 _sameFilter    ( false),
  95 _rejectPileup  ( 1),
  96 _rejectPairConversion ( 0),
  97 _vertexZMin           ( -10),
  98 _vertexZMax           (  10),
  99 _vertexXYMin          ( -10),
 100 _vertexXYMax          (  10),
 101 _centralityMethod     (  4),
 102 _centralityMin        (  0.),
 103 _centralityMax        (  0.),
 104 _requestedCharge_1    (   1),
 105 _requestedCharge_2    (  -1),
 106 _dcaZMin              ( -3),
 107 _dcaZMax              (  3.),
 108 _dcaXYMin             ( -2.4),
 109 _dcaXYMax             (  2.4),
 110 _dedxMin              ( 0),
 111 _dedxMax              ( 100000),
 112 _nClusterMin          ( 80),
 113 _trackFilterBit       (0),
 114 fNSigmaCut            (3.),
 115 //_min_eta_1(0),
 116 //_max_eta_1 (0),
 117 //_min_eta_2 (0),
 118 //_max_eta_2 (0),
 119
 120 _field    ( 1.),
 121 _nTracks  ( 0 ),
 122 _mult0    ( 0 ),
 123 _mult1    ( 0 ),
 124 _mult2    ( 0 ),
 125 _mult3    ( 0 ),
 126 _mult4    ( 0 ),
 127 _mult4a    ( 0 ),
 128 _mult5    ( 0 ),
 129 _mult6    ( 0 ),
 130 arraySize ( 2000),
 131 _id_1(0),
 132 _charge_1(0),
 133 _iEtaPhi_1(0),
 134 _iPt_1(0),
 135 _pt_1(0),
 136 _px_1(0),
 137 _py_1(0),
 138 _pz_1(0),
 139 _correction_1(0),
 140 _dedx_1(0),
 141 _id_2(0),
 142 _charge_2(0),
 143 _iEtaPhi_2(0),
 144 _iPt_2(0),
 145 _pt_2(0),
 146 _px_2(0),
 147 _py_2(0),
 148 _pz_2(0),
 149 _correction_2(0),
 150 _dedx_2(0),
 151 _correctionWeight_1(0),
 152 _correctionWeight_2(0),
 153 _nBins_M0(500),       _min_M0(0),        _max_M0(10000),          _width_M0(20),
 154 _nBins_M1(500),       _min_M1(0),        _max_M1(10000),          _width_M1(20),
 155 _nBins_M2(500),       _min_M2(0),        _max_M2(10000),          _width_M2(20),
 156 _nBins_M3(500),       _min_M3(0),        _max_M3(10000),          _width_M3(20),
 157 _nBins_M4(100),       _min_M4(0),        _max_M4(1),              _width_M4(0.01),
 158 _nBins_M5(100),       _min_M5(0),        _max_M5(1),              _width_M5(0.01),
 159 _nBins_M6(100),       _min_M6(0),        _max_M6(1),              _width_M6(0.01),
 160 _nBins_vertexZ(40),   _min_vertexZ(-10), _max_vertexZ(10),        _width_vertexZ(0.5),
 161
 162 _nBins_pt_1(18),      _min_pt_1(0.2),    _max_pt_1(2.0),          _width_pt_1(0.1),
 163 _nBins_phi_1(72),     _min_phi_1(0),     _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
 164 //_nBins_eta_1(18),
 165 //_min_eta_1(-0.9),  _max_eta_1(0.9),           _width_eta_1(0.1),
 166   _nBins_eta_1(0),
 167   _min_eta_1(0),  _max_eta_1(0),           _width_eta_1(0.1),
 168
 169 _nBins_etaPhi_1(0),
 170 _nBins_etaPhiPt_1(0),
 171 _nBins_zEtaPhiPt_1(0),
 172
 173 _nBins_pt_2(18),     _min_pt_2(0.2),     _max_pt_2(2.0),          _width_pt_2(0.1),
 174 _nBins_phi_2(72),    _min_phi_2(0),      _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
 175 //_nBins_eta_2(18),
 176 //_min_eta_2(-0.9),     _max_eta_2(0.9),           _width_eta_2(0.1),
 177   _nBins_eta_2(0),
 178   _min_eta_2(0),     _max_eta_2(0),           _width_eta_2(0.1),
 179
 180
 181
 182 _nBins_etaPhi_2(0),
 183 _nBins_etaPhiPt_2(0),
 184 _nBins_zEtaPhiPt_2(0),
 185 _nBins_etaPhi_12(0),
 186 __n1_1(0),
 187 __n1_2(0),
 188 __n2_12(0),
 189 __s1pt_1(0),
 190 __s1pt_2(0),
 191 __s2ptpt_12(0),
 192 __s2NPt_12(0),
 193 __s2PtN_12(0),
 194 __n1Nw_1(0),
 195 __n1Nw_2(0),
 196 __n2Nw_12(0),
 197 __s1ptNw_1(0),
 198 __s1ptNw_2(0),
 199 __s2ptptNw_12(0),
 200 __s2NPtNw_12(0),
 201 __s2PtNNw_12(0),
 202 __n1_1_vsPt(0),
 203 __n1_1_vsEtaPhi(0),
 204 __s1pt_1_vsEtaPhi(0),
 205 __n1_1_vsZEtaPhiPt(0),
 206 __n1_2_vsPt(0),
 207 __n1_2_vsEtaPhi(0),
 208 __s1pt_2_vsEtaPhi(0),
 209 __n1_2_vsZEtaPhiPt(0),
 210 __n2_12_vsPtPt(0),
 211 __n2_12_vsEtaPhi(0),
 212 __s2ptpt_12_vsEtaPhi(0),
 213 __s2PtN_12_vsEtaPhi(0),
 214 __s2NPt_12_vsEtaPhi(0),
 215 _weight_1      ( 0    ),
 216 _weight_2      ( 0    ),
 217 _eventAccounting ( 0),
 218 _m0 ( 0),
 219 _m1 ( 0),
 220 _m2 ( 0),
 221 _m3 ( 0),
 222 _m4 ( 0),
 223 _m5 ( 0),
 224 _m6 ( 0),
 225 _vertexZ ( 0),
 226
 227 _Ncluster1  ( 0),
 228 _Ncluster2  ( 0),
 229 _etadis ( 0),
 230 _phidis ( 0),
 231 _dcaz   ( 0),
 232 _dcaxy  ( 0),
 233 _n1_1_vsPt         ( 0),
 234 _n1_1_vsEtaVsPhi   ( 0),
 235 _s1pt_1_vsEtaVsPhi ( 0),
 236 _n1_1_vsZVsEtaVsPhiVsPt ( 0),
 237 _n1_1_vsM          ( 0),  // w/ weight
 238 _s1pt_1_vsM        ( 0),
 239 _n1Nw_1_vsM        ( 0), // w/o weight
 240 _s1ptNw_1_vsM      ( 0),
 241 _dedxVsP_1         ( 0),
 242 _corrDedxVsP_1     ( 0),
 243 _betaVsP_1         ( 0),
 244 _n1_2_vsPt         ( 0),
 245 _n1_2_vsEtaVsPhi   ( 0),
 246 _s1pt_2_vsEtaVsPhi ( 0),
 247 _n1_2_vsZVsEtaVsPhiVsPt ( 0),
 248 _n1_2_vsM          ( 0),
 249 _s1pt_2_vsM        ( 0),
 250 _n1Nw_2_vsM        ( 0),
 251 _s1ptNw_2_vsM      ( 0),
 252 _dedxVsP_2         ( 0),
 253 _corrDedxVsP_2     ( 0),
 254 _betaVsP_2         ( 0),
 255 _n2_12_vsEtaPhi    ( 0),
 256 _n2_12_vsPtVsPt    ( 0),
 257 _s2PtPt_12_vsEtaPhi( 0),
 258 _s2PtN_12_vsEtaPhi ( 0),
 259 _s2NPt_12_vsEtaPhi ( 0),
 260 _n2_12_vsM         ( 0),
 261 _s2PtPt_12_vsM     ( 0),
 262 _s2PtN_12_vsM      ( 0),
 263 _s2NPt_12_vsM      ( 0),
 264 _n2Nw_12_vsM       ( 0),
 265 _s2PtPtNw_12_vsM   ( 0),
 266 _s2PtNNw_12_vsM    ( 0),
 267 _s2NPtNw_12_vsM    ( 0),
 268 _invMass           ( 0),
 269 _invMassElec       ( 0),
 270 n1Name("NA"),
 271 n1NwName("NA"),
 272 n2Name("NA"),
 273 n2NwName("NA"),
 274 n3Name("NA"),
 275 n1n1Name("NA"),
 276 n1n1n1Name("NA"),
 277 n2n1Name("NA"),
 278 r1Name("NA"),
 279 r2Name("NA"),
 280 r3Name("NA"),
 281 r2r1Name("NA"),
 282 c2Name("NA"),
 283 c3Name("NA"),
 284 d3Name("NA"),
 285 p3Name("NA"),
 286 cName("NA"),
 287
 288 intR2Name("NA"),
 289 binCorrName("NA"),
 290 intBinCorrName("NA"),
 291
 292 countsName("NA"),
 293 part_1_Name("NA"),
 294 part_2_Name("NA"),
 295 part_3_Name("NA"),
 296 pair_12_Name("NA"),
 297 pair_13_Name("NA"),
 298 pair_23_Name("NA"),
 299 tripletName("NA"),
 300
 301 avg("NA"),
 302 avgName("NA"),
 303 sumName("NA"),
 304 s1ptName("NA"),
 305 s1ptNwName("NA"),
 306 s1DptName("NA"),
 307
 308 s2PtPtName("NA"),
 309 s2NPtName("NA"),
 310 s2PtNName("NA"),
 311 s2DptDptName("NA"),
 312
 313 s2PtPtNwName("NA"),
 314 s2NPtNwName("NA"),
 315 s2PtNNwName("NA"),
 316
 317 ptName("NA"),
 318 ptptName("NA"),
 319 pt1pt1Name("NA"),
 320 DptName("NA"),
 321 DptDptName("NA"),
 322 RDptDptName("NA"),
 323 nPtName("NA"),
 324 ptNName("NA"),
 325 seanName("NA"),
 326
 327 _title_counts("NA"),
 328
 329 _title_m0("NA"),
 330 _title_m1("NA"),
 331 _title_m2("NA"),
 332 _title_m3("NA"),
 333 _title_m4("NA"),
 334 _title_m5("NA"),
 335 _title_m6("NA"),
 336
 337 _title_eta_1("NA"),
 338 _title_phi_1("NA"),
 339 _title_pt_1("NA"),
 340 _title_etaPhi_1("NA"),
 341 _title_n_1("NA"),
 342 _title_SumPt_1("NA"),
 343 _title_AvgPt_1("NA"),
 344 _title_AvgN_1("NA"),
 345 _title_AvgSumPt_1("NA"),
 346
 347 _title_eta_2("NA"),
 348 _title_phi_2("NA"),
 349 _title_pt_2("NA"),
 350 _title_etaPhi_2("NA"),
 351 _title_n_2("NA"),
 352 _title_SumPt_2("NA"),
 353 _title_AvgPt_2("NA"),
 354 _title_AvgN_2("NA"),
 355 _title_AvgSumPt_2("NA"),
 356
 357 _title_etaPhi_12("NA"),
 358
 359 _title_AvgN2_12("NA"),
 360 _title_AvgSumPtPt_12("NA"),
 361 _title_AvgSumPtN_12("NA"),
 362 _title_AvgNSumPt_12("NA"),
 363
 364 vsZ("NA"),
 365 vsM("NA"),
 366 vsPt("NA"),
 367 vsPhi("NA"),
 368 vsEta("NA"),
 369 vsEtaPhi("NA"),
 370 vsPtVsPt("NA")
 371 {
 372   printf("Default constructor called \n");
 373
 374   printf("passed \n ");
 375
 376 }
 377
 378 AliAnalysisTaskDptDptCorrelations::AliAnalysisTaskDptDptCorrelations(const TString & name)
 379 : AliAnalysisTaskSE(name),
 380 fAODEvent(0),
 381 fESDEvent(0),
 382 fInputHandler(0),
 383 fPIDResponse(0),
 384 _outputHistoList(0),
 385 _twoPi         ( 2.0 * 3.1415927),
 386 _eventCount    ( 0),
 387 _debugLevel    ( 0),
 388 _singlesOnly   ( 0),
 389 _useWeights    ( 0),
 390 _sameFilter    ( false),
 391 _rejectPileup  ( 1),
 392 _rejectPairConversion ( 0),
 393 _vertexZMin           ( -10.),
 394 _vertexZMax           (  10.),
 395 _vertexXYMin          ( -10.),
 396 _vertexXYMax          (  10.),
 397 _centralityMethod     (  4),
 398 _centralityMin        (  0.),
 399 _centralityMax        (  1.),
 400 _requestedCharge_1    (   1),
 401 _requestedCharge_2    (  -1),
 402 _dcaZMin              ( -3),
 403 _dcaZMax              (  3.),
 404 _dcaXYMin             ( -2.4),
 405 _dcaXYMax             (  2.4),
 406 _dedxMin              ( 0),
 407 _dedxMax              ( 100000),
 408 _nClusterMin          ( 80),
 409 _trackFilterBit       ( 0),
 410 fNSigmaCut            ( 3.),
 411 _field    ( 1.),
 412 _nTracks  ( 0 ),
 413 _mult0    ( 0 ),
 414 _mult1    ( 0 ),
 415 _mult2    ( 0 ),
 416 _mult3    ( 0 ),
 417 _mult4    ( 0 ),
 418 _mult4a    ( 0 ),
 419 _mult5    ( 0 ),
 420 _mult6    ( 0 ),
 421 arraySize ( 2000),
 422 _id_1(0),
 423 _charge_1(0),
 424 _iEtaPhi_1(0),
 425 _iPt_1(0),
 426 _pt_1(0),
 427 _px_1(0),
 428 _py_1(0),
 429 _pz_1(0),
 430 _correction_1(0),
 431 _dedx_1(0),
 432 _id_2(0),
 433 _charge_2(0),
 434 _iEtaPhi_2(0),
 435 _iPt_2(0),
 436 _pt_2(0),
 437 _px_2(0),
 438 _py_2(0),
 439 _pz_2(0),
 440 _correction_2(0),
 441 _dedx_2(0),
 442 _correctionWeight_1(0),
 443 _correctionWeight_2(0),
 444 _nBins_M0(500),       _min_M0(0),        _max_M0(10000),          _width_M0(20),
 445 _nBins_M1(500),       _min_M1(0),        _max_M1(10000),          _width_M1(20),
 446 _nBins_M2(500),       _min_M2(0),        _max_M2(10000),          _width_M2(20),
 447 _nBins_M3(500),       _min_M3(0),        _max_M3(10000),          _width_M3(20),
 448 _nBins_M4(100),       _min_M4(0),        _max_M4(1),              _width_M4(0.01),
 449 _nBins_M5(100),       _min_M5(0),        _max_M5(1),              _width_M5(0.01),
 450 _nBins_M6(100),       _min_M6(0),        _max_M6(1),              _width_M6(0.01),
 451 _nBins_vertexZ(40),   _min_vertexZ(-10), _max_vertexZ(10),        _width_vertexZ(0.5),
 452
 453 _nBins_pt_1(18),      _min_pt_1(0.2),    _max_pt_1(2.0),          _width_pt_1(0.1),
 454 _nBins_phi_1(72),     _min_phi_1(0),     _max_phi_1(2.*3.1415927),_width_phi_1(2.*3.1415927/72.),
 455   //_nBins_eta_1(18),
 456   //_min_eta_1(-0.9),    _max_eta_1(0.9),           _width_eta_1(0.1),
 457   _nBins_eta_1(0),
 458   _min_eta_1(0),    _max_eta_1(0),           _width_eta_1(0.1),
 459
 460 _nBins_etaPhi_1(0),
 461 _nBins_etaPhiPt_1(0),
 462 _nBins_zEtaPhiPt_1(0),
 463
 464 _nBins_pt_2(18),     _min_pt_2(0.2),     _max_pt_2(2.0),          _width_pt_2(0.1),
 465 _nBins_phi_2(72),    _min_phi_2(0),      _max_phi_2(2.*3.1415927),_width_phi_2(2.*3.1415927/72),
 466   //_nBins_eta_2(18),
 467   //_min_eta_2(-0.9),     _max_eta_2(0.9),           _width_eta_2(0.1),
 468   _nBins_eta_2(0),
 469   _min_eta_2(0),     _max_eta_2(0),           _width_eta_2(0.1),
 470
 471
 472 _nBins_etaPhi_2(0),
 473 _nBins_etaPhiPt_2(0),
 474 _nBins_zEtaPhiPt_2(0),
 475 _nBins_etaPhi_12(0),
 476 __n1_1(0),
 477 __n1_2(0),
 478 __n2_12(0),
 479 __s1pt_1(0),
 480 __s1pt_2(0),
 481 __s2ptpt_12(0),
 482 __s2NPt_12(0),
 483 __s2PtN_12(0),
 484 __n1Nw_1(0),
 485 __n1Nw_2(0),
 486 __n2Nw_12(0),
 487 __s1ptNw_1(0),
 488 __s1ptNw_2(0),
 489 __s2ptptNw_12(0),
 490 __s2NPtNw_12(0),
 491 __s2PtNNw_12(0),
 492 __n1_1_vsPt(0),
 493 __n1_1_vsEtaPhi(0),
 494 __s1pt_1_vsEtaPhi(0),
 495 __n1_1_vsZEtaPhiPt(0),
 496 __n1_2_vsPt(0),
 497 __n1_2_vsEtaPhi(0),
 498 __s1pt_2_vsEtaPhi(0),
 499 __n1_2_vsZEtaPhiPt(0),
 500 __n2_12_vsPtPt(0),
 501 __n2_12_vsEtaPhi(0),
 502 __s2ptpt_12_vsEtaPhi(0),
 503 __s2PtN_12_vsEtaPhi(0),
 504 __s2NPt_12_vsEtaPhi(0),
 505 _weight_1        ( 0    ),
 506 _weight_2        ( 0    ),
 507 _eventAccounting ( 0),
 508 _m0 ( 0),
 509 _m1 ( 0),
 510 _m2 ( 0),
 511 _m3 ( 0),
 512 _m4 ( 0),
 513 _m5 ( 0),
 514 _m6 ( 0),
 515 _vertexZ ( 0),
 516 _Ncluster1  ( 0),
 517 _Ncluster2  ( 0),
 518 _etadis ( 0),
 519 _phidis ( 0),
 520
 521 _dcaz ( 0),
 522 _dcaxy ( 0),
 523 _n1_1_vsPt         ( 0),
 524 _n1_1_vsEtaVsPhi   ( 0),
 525 _s1pt_1_vsEtaVsPhi ( 0),
 526 _n1_1_vsZVsEtaVsPhiVsPt ( 0),
 527 _n1_1_vsM          ( 0),  // w/ weight
 528 _s1pt_1_vsM        ( 0),
 529 _n1Nw_1_vsM        ( 0), // w/o weight
 530 _s1ptNw_1_vsM      ( 0),
 531 _dedxVsP_1         ( 0),
 532 _corrDedxVsP_1     ( 0),
 533 _betaVsP_1         ( 0),
 534 _n1_2_vsPt         ( 0),
 535 _n1_2_vsEtaVsPhi   ( 0),
 536 _s1pt_2_vsEtaVsPhi ( 0),
 537 _n1_2_vsZVsEtaVsPhiVsPt ( 0),
 538 _n1_2_vsM          ( 0),
 539 _s1pt_2_vsM        ( 0),
 540 _n1Nw_2_vsM        ( 0),
 541 _s1ptNw_2_vsM      ( 0),
 542 _dedxVsP_2         ( 0),
 543 _corrDedxVsP_2     ( 0),
 544 _betaVsP_2         ( 0),
 545 _n2_12_vsEtaPhi    ( 0),
 546 _n2_12_vsPtVsPt    ( 0),
 547 _s2PtPt_12_vsEtaPhi( 0),
 548 _s2PtN_12_vsEtaPhi ( 0),
 549 _s2NPt_12_vsEtaPhi ( 0),
 550 _n2_12_vsM         ( 0),
 551 _s2PtPt_12_vsM     ( 0),
 552 _s2PtN_12_vsM      ( 0),
 553 _s2NPt_12_vsM      ( 0),
 554 _n2Nw_12_vsM       ( 0),
 555 _s2PtPtNw_12_vsM   ( 0),
 556 _s2PtNNw_12_vsM    ( 0),
 557 _s2NPtNw_12_vsM    ( 0),
 558 _invMass           ( 0),
 559 _invMassElec       ( 0),
 560 n1Name("NA"),
 561 n1NwName("NA"),
 562 n2Name("NA"),
 563 n2NwName("NA"),
 564 n3Name("NA"),
 565 n1n1Name("NA"),
 566 n1n1n1Name("NA"),
 567 n2n1Name("NA"),
 568 r1Name("NA"),
 569 r2Name("NA"),
 570 r3Name("NA"),
 571 r2r1Name("NA"),
 572 c2Name("NA"),
 573 c3Name("NA"),
 574 d3Name("NA"),
 575 p3Name("NA"),
 576 cName("NA"),
 577
 578 intR2Name("NA"),
 579 binCorrName("NA"),
 580 intBinCorrName("NA"),
 581
 582 countsName("NA"),
 583 part_1_Name("NA"),
 584 part_2_Name("NA"),
 585 part_3_Name("NA"),
 586 pair_12_Name("NA"),
 587 pair_13_Name("NA"),
 588 pair_23_Name("NA"),
 589 tripletName("NA"),
 590
 591 avg("NA"),
 592 avgName("NA"),
 593 sumName("NA"),
 594 s1ptName("NA"),
 595 s1ptNwName("NA"),
 596 s1DptName("NA"),
 597
 598 s2PtPtName("NA"),
 599 s2NPtName("NA"),
 600 s2PtNName("NA"),
 601 s2DptDptName("NA"),
 602
 603 s2PtPtNwName("NA"),
 604 s2NPtNwName("NA"),
 605 s2PtNNwName("NA"),
 606
 607 ptName("NA"),
 608 ptptName("NA"),
 609 pt1pt1Name("NA"),
 610 DptName("NA"),
 611 DptDptName("NA"),
 612 RDptDptName("NA"),
 613 nPtName("NA"),
 614 ptNName("NA"),
 615 seanName("NA"),
 616
 617 _title_counts("NA"),
 618
 619 _title_m0("NA"),
 620 _title_m1("NA"),
 621 _title_m2("NA"),
 622 _title_m3("NA"),
 623 _title_m4("NA"),
 624 _title_m5("NA"),
 625 _title_m6("NA"),
 626
 627 _title_eta_1("NA"),
 628 _title_phi_1("NA"),
 629 _title_pt_1("NA"),
 630 _title_etaPhi_1("NA"),
 631 _title_n_1("NA"),
 632 _title_SumPt_1("NA"),
 633 _title_AvgPt_1("NA"),
 634 _title_AvgN_1("NA"),
 635 _title_AvgSumPt_1("NA"),
 636
 637 _title_eta_2("NA"),
 638 _title_phi_2("NA"),
 639 _title_pt_2("NA"),
 640 _title_etaPhi_2("NA"),
 641 _title_n_2("NA"),
 642 _title_SumPt_2("NA"),
 643 _title_AvgPt_2("NA"),
 644 _title_AvgN_2("NA"),
 645 _title_AvgSumPt_2("NA"),
 646
 647 _title_etaPhi_12("NA"),
 648
 649 _title_AvgN2_12("NA"),
 650 _title_AvgSumPtPt_12("NA"),
 651 _title_AvgSumPtN_12("NA"),
 652 _title_AvgNSumPt_12("NA"),
 653
 654 vsZ("NA"),
 655 vsM("NA"),
 656 vsPt("NA"),
 657 vsPhi("NA"),
 658 vsEta("NA"),
 659 vsEtaPhi("NA"),
 660 vsPtVsPt("NA")
 661 {
 662   printf("2nd constructor called ");
 663
 664   DefineOutput(0, TList::Class());
 665
 666   printf("passed  ");
 667
 668 }
 669
 670 AliAnalysisTaskDptDptCorrelations::~AliAnalysisTaskDptDptCorrelations()
 671 {
 672   /*
 673    delete _id_1;
 674    delete _charge_1;
 675    delete _iEtaPhi_1;
 676    delete _iPt_1;
 677    delete _pt_1;
 678    delete _px_1;
 679    delete _py_1;
 680    delete _pz_1;
 681    delete _correction_1;
 682    delete _dedx_1;
 683    delete __n1_1_vsPt;
 684    delete __n1_1_vsEtaPhi;
 685    delete __s1pt_1_vsEtaPhi;
 686    delete __n1_1_vsZEtaPhiPt;
 687    if (_correctionWeight_1) delete _correctionWeight_1;
 688
 689    if (!_sameFilter)
 690    {
 691    delete _id_2;
 692    delete _charge_2;
 693    delete _iEtaPhi_2;
 694    delete _iPt_2;
 695    delete _pt_2;
 696    delete _px_2;
 697    delete _py_2;
 698    delete _pz_2;
 699    delete _correction_2;
 700    delete _dedx_2;
 701    delete __n1_2_vsPt;
 702    delete __n1_2_vsEtaPhi;
 703    delete __s1pt_2_vsEtaPhi;
 704    delete __n1_2_vsZEtaPhiPt;
 705    if (_correctionWeight_2) delete _correctionWeight_2;
 706    }
 707
 708    if (!_singlesOnly)
 709    {
 710    delete __n2_12_vsPtPt;
 711    delete __n2_12_vsEtaPhi;
 712    delete __s2ptpt_12_vsEtaPhi;
 713    delete __s2PtN_12_vsEtaPhi;
 714    delete __s2NPt_12_vsEtaPhi;
 715    }
 716    */
 717 }
 718
 719 void AliAnalysisTaskDptDptCorrelations::UserCreateOutputObjects()
 720 {
 721   //cout<< "AliAnalysisTaskDptDptCorrelations::CreateOutputObjects() Starting " << endl;
 722   OpenFile(0);
 723   _outputHistoList = new TList();
 724   _outputHistoList->SetOwner();
 725
 726   //if (_useWeights) DefineInput(2, TList::Class());
 727
 728   //Setup the parameters of the histograms
 729   _nBins_M0 = 500; _min_M0   = 0.;    _max_M0    = 5000.;  _width_M0 = (_max_M0-_min_M0)/_nBins_M0;
 730   _nBins_M1 = 500; _min_M1   = 0.;    _max_M1    = 5000.;  _width_M1 = (_max_M1-_min_M1)/_nBins_M1;
 731   _nBins_M2 = 500; _min_M2   = 0.;    _max_M2    = 5000.;  _width_M2 = (_max_M2-_min_M2)/_nBins_M2;
 732   _nBins_M3 = 500; _min_M3   = 0.;    _max_M3    = 5000.;  _width_M3 = (_max_M3-_min_M3)/_nBins_M3;
 733   _nBins_M4 = 100; _min_M4   = 0.;    _max_M4    = 100.;   _width_M4 = (_max_M4-_min_M4)/_nBins_M4;
 734   _nBins_M5 = 100; _min_M5   = 0.;    _max_M5    = 100.;   _width_M5 = (_max_M5-_min_M5)/_nBins_M5;
 735   _nBins_M6 = 100; _min_M6   = 0.;    _max_M6    = 100.;   _width_M6 = (_max_M6-_min_M6)/_nBins_M6;
 736
 737   _min_vertexZ       = _vertexZMin;
 738   _max_vertexZ       = _vertexZMax;
 739   _width_vertexZ     = 0.5;
 740   _nBins_vertexZ     = int(0.5+ (_max_vertexZ - _min_vertexZ)/_width_vertexZ);
 741   _nBins_pt_1        = int(0.5+ (_max_pt_1 -_min_pt_1 )/_width_pt_1);
 742   _nBins_eta_1       = int(0.5+ (_max_eta_1-_min_eta_1)/_width_eta_1);
 743   _width_phi_1       = (_max_phi_1  - _min_phi_1)  /_nBins_phi_1;
 744   _nBins_etaPhi_1    = _nBins_phi_1    * _nBins_eta_1;
 745   _nBins_etaPhiPt_1  = _nBins_etaPhi_1 * _nBins_pt_1;
 746   _nBins_zEtaPhiPt_1 = _nBins_vertexZ  * _nBins_etaPhiPt_1;
 747
 748   _nBins_pt_2        =  int(0.5+ (_max_pt_2 -_min_pt_2 )/_width_pt_2);
 749   _nBins_eta_2       = int(0.5+ (_max_eta_2-_min_eta_2)/_width_eta_2);
 750   _width_phi_2       = (_max_phi_2  - _min_phi_2)  /_nBins_phi_2;
 751   _nBins_etaPhi_2    = _nBins_phi_2    * _nBins_eta_2;
 752   _nBins_etaPhiPt_2  = _nBins_etaPhi_2 * _nBins_pt_2;
 753   _nBins_zEtaPhiPt_2 = _nBins_vertexZ  * _nBins_etaPhiPt_2;
 754   _nBins_etaPhi_12   = _nBins_etaPhi_1 * _nBins_etaPhi_2;
 755
 756   //setup the work arrays
 757
 758   _id_1       = new int[arraySize];
 759   _charge_1   = new int[arraySize];
 760   //_iPhi_1 = new int[arraySize];
 761   //_iEta_1 = new int[arraySize];
 762   _iEtaPhi_1  = new int[arraySize];
 763   _iPt_1      = new int[arraySize];
 764   _pt_1       = new float[arraySize];
 765   _px_1       = new float[arraySize];
 766   _py_1       = new float[arraySize];
 767   _pz_1       = new float[arraySize];
 768   //_phi_1 = new float[arraySize];
 769   //_eta_1 = new float[arraySize];
 770   _correction_1 = new float[arraySize];
 771   //_dedx_1     = new float[arraySize];
 772
 773   __n1_1_vsPt              = getDoubleArray(_nBins_pt_1,        0.);
 774   __n1_1_vsEtaPhi          = getDoubleArray(_nBins_etaPhi_1,    0.);
 775   __s1pt_1_vsEtaPhi        = getDoubleArray(_nBins_etaPhi_1,    0.);
 776   __n1_1_vsZEtaPhiPt       = getFloatArray(_nBins_zEtaPhiPt_1,  0.);
 777
 778   //cout << "==========================================================================================" << endl;
 779   //cout << "=============== Booking for particle 1 done." << endl;
 780   //cout << "_requestedCharge_1: " << _requestedCharge_1 << endl;
 781   //cout << "_requestedCharge_2: " << _requestedCharge_2 << endl;
 782
 783   if (_requestedCharge_2!=_requestedCharge_1)
 784     {
 785       //cout << " creating arrays for particle 2 with size: " << arraySize << endl;
 786     _sameFilter = 0;
 787     //particle 2
 788     _id_2       = new int[arraySize];
 789     _charge_2   = new int[arraySize];
 790     //_iPhi_2   = new int[arraySize];
 791     //_iEta_2   = new int[arraySize];
 792     _iEtaPhi_2  = new int[arraySize];
 793     _iPt_2      = new int[arraySize];
 794     _pt_2       = new float[arraySize];
 795     _px_2       = new float[arraySize];
 796     _py_2       = new float[arraySize];
 797     _pz_2       = new float[arraySize];
 798     //_phi_2    = new float[arraySize];
 799     //_eta_2    = new float[arraySize];
 800     _correction_2 = new float[arraySize];
 801     //_dedx_2       = new float[arraySize];
 802
 803     __n1_2_vsPt              = getDoubleArray(_nBins_pt_2,        0.);
 804     __n1_2_vsEtaPhi          = getDoubleArray(_nBins_etaPhi_2,    0.);
 805     __s1pt_2_vsEtaPhi        = getDoubleArray(_nBins_etaPhi_2,    0.);
 806     __n1_2_vsZEtaPhiPt       = getFloatArray(_nBins_zEtaPhiPt_2, 0.);
 807
 808     }
 809
 810   __n2_12_vsPtPt           = getDoubleArray(_nBins_pt_1*_nBins_pt_2,0.);
 811   __n2_12_vsEtaPhi         = getFloatArray(_nBins_etaPhi_12,       0.);
 812   __s2ptpt_12_vsEtaPhi     = getFloatArray(_nBins_etaPhi_12,       0.);
 813   __s2PtN_12_vsEtaPhi      = getFloatArray(_nBins_etaPhi_12,       0.);
 814   __s2NPt_12_vsEtaPhi      = getFloatArray(_nBins_etaPhi_12,       0.);
 815
 816   // Setup all the labels needed.
 817
 818   part_1_Name   = "_1";
 819   part_2_Name   = "_2";
 820   pair_12_Name  = "_12";
 821
 822   n1Name     = "n1";
 823   n2Name     = "n2";
 824   n1NwName   = "n1Nw";
 825   n2NwName   = "n2Nw";
 826   r1Name     = "r1";
 827   r2Name     = "r2";
 828   r3Name     = "r3";
 829   r2r1Name   = "r2r1";
 830   c2Name     = "c2";
 831   c3Name     = "c3";
 832   d3Name     = "d3";
 833   p3Name     = "p3";
 834   cName      = "sean";
 835
 836   intR2Name       = "intR2";
 837   binCorrName     = "binCorr";
 838   intBinCorrName  = "intBinCorr";
 839
 840   avgName      = "avg";
 841   sumName      = "sum";
 842   s1ptName     = "sumPt";
 843   s1ptNwName   = "sumPtNw";
 844   s1DptName    = "sumDpt";
 845   s2PtPtName   = "sumPtPt";
 846   s2PtPtNwName = "sumPtPtNw";
 847   s2DptDptName = "sumDptDpt";
 848   s2NPtName    = "sumNPt";
 849   s2NPtNwName  = "sumNPtNw";
 850   s2PtNName    = "sumPtN";
 851   s2NPtNwName  = "sumNPtNw";
 852   s2PtNNwName  = "sumPtNNw";
 853   ptName       = "avgPt";
 854   ptptName     = "avgPtPt";
 855   pt1pt1Name   = "avgPtavgPt";
 856   DptName      = "avgDpt";
 857   DptDptName   = "avgDptDpt";
 858   RDptDptName  = "relDptDpt"; // ratio of avgDptDpt by avgPt*avgPt
 859   nPtName      = "avgNpt";
 860   ptNName      = "avgPtN";
 861   seanName     = "seanC";
 862
 863   _title_counts = "yield";
 864
 865   _title_m0     = "M_{0}";
 866   _title_m1     = "M_{1}";
 867   _title_m2     = "M_{2}";
 868   _title_m3     = "M_{3}";
 869   _title_m4     = "V0Centrality";
 870   _title_m5     = "TrkCentrality";
 871   _title_m6     = "SpdCentrality";
 872
 873   _title_eta_1       = "#eta_{1}";
 874   _title_phi_1       = "#varphi_{1} (radian)";
 875   _title_etaPhi_1    = "#eta_{1}#times#varphi_{1}";
 876   _title_pt_1        = "p_{t,1} (GeV/c)";
 877   _title_n_1         = "n_{1}";
 878   _title_SumPt_1     = "#Sigma p_{t,1} (GeV/c)";
 879   _title_AvgPt_1     = "#LT p_{t,1} #GT (GeV/c)";
 880   _title_AvgN_1      = "#LT n_{1} #GT";
 881   _title_AvgSumPt_1  = "#LT #Sigma p_{t,1} #GT (GeV/c)";
 882
 883   _title_eta_2       = "#eta_{2}";
 884   _title_phi_2       = "#varphi_{2} (radian)";
 885   _title_etaPhi_2    = "#eta_{2}#times#varphi_{2}";
 886   _title_pt_2        = "p_{t,2} (GeV/c)";
 887   _title_n_2         = "n_{2}";
 888   _title_SumPt_2     = "#Sigma p_{t,1} (GeV/c)";
 889   _title_AvgPt_2     = "#LT p_{t,2} #GT (GeV/c)";
 890   _title_AvgN_2      = "#LT n_{2} #GT";
 891   _title_AvgSumPt_2  = "#LT #Sigma p_{t,2} #GT (GeV/c)";
 892
 893   _title_etaPhi_12   = "#eta_{1}#times#varphi_{1}#times#eta_{2}#times#varphi_{2}";
 894
 895   _title_AvgN2_12       = "#LT n_{2} #GT";;
 896   _title_AvgSumPtPt_12  = "#LT #Sigma p_{t,1}p_{t,2} #GT";;
 897   _title_AvgSumPtN_12   = "#LT #Sigma p_{t,1}N #GT";;
 898   _title_AvgNSumPt_12   = "#LT N#Sigma p_{t,2} #GT";;
 899
 900
 901   vsZ         = "_vsZ";
 902   vsM         = "_vsM";
 903   vsPt         = "_vsPt";
 904   vsPhi        = "_vsPhi";
 905   vsEta        = "_vsEta";
 906   vsEtaPhi     = "_vsEtaPhi";
 907   vsPtVsPt     = "_vsPtVsPt";
 908
 909
 910   if (_useWeights)
 911     {
 912     int iZ, iEtaPhi, iPt;
 913     int iZ1,iEtaPhi1,iPt1;
 914     int a, b;
 915     if (_weight_1)
 916       {
 917       _correctionWeight_1 = new float[_nBins_vertexZ*_nBins_etaPhi_1*_nBins_pt_1];
 918       a = _nBins_pt_1;
 919       b = _nBins_etaPhi_1*_nBins_pt_1;
 920       for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
 921         {
 922         for (iEtaPhi=0,iEtaPhi1=1; iEtaPhi<_nBins_etaPhi_1; iEtaPhi++, iEtaPhi1++)
 923           {
 924           for (iPt=0,iPt1=1; iPt<_nBins_pt_1; iPt++, iPt1++)
 925             {
 926             _correctionWeight_1[iZ*b+iEtaPhi*a+iPt] = _weight_1->GetBinContent(iZ1,iEtaPhi1,iPt1);
 927             }
 928           }
 929         }
 930       } // _weight_1
 931     else
 932       {
 933       AliError("AliAnalysisTaskDptDptCorrelations:: _weight_1 is a null pointer.");
 934       return;
 935       }
 936     if (!_sameFilter)
 937       {
 938       if (_weight_2)
 939         {
 940         _correctionWeight_2 = new float[_nBins_vertexZ*_nBins_etaPhi_2*_nBins_pt_2];
 941         a = _nBins_pt_2;
 942         b = _nBins_etaPhi_2*_nBins_pt_2;
 943         for (iZ=0,iZ1=1; iZ<_nBins_vertexZ; iZ++, iZ1++)
 944           {
 945           for (iEtaPhi=0,iEtaPhi1=1; iEtaPhi<_nBins_etaPhi_2; iEtaPhi++, iEtaPhi1++)
 946             {
 947             for (iPt=0,iPt1=1; iPt<_nBins_pt_2; iPt++, iPt1++)
 948               {
 949               _correctionWeight_2[iZ*b+iEtaPhi*a+iPt] = _weight_2->GetBinContent(iZ1,iEtaPhi1,iPt1);
 950               }
 951             }
 952           }
 953         } // _weight_2
 954       else
 955         {
 956         AliError("AliAnalysisTaskDptDptCorrelations:: _weight_1 is a null pointer.");
 957         return;
 958         }
 959       }
 960     }
 961
 962   createHistograms();
 963   PostData(0,_outputHistoList);
 964
 965   //cout<< "AliAnalysisTaskDptDptCorrelations::CreateOutputObjects() DONE " << endl;
 966
 967 }
 968
 969 void  AliAnalysisTaskDptDptCorrelations::createHistograms()
 970 {
 971   AliInfo(" AliAnalysisTaskDptDptCorrelations::createHistoHistograms() Creating Event Histos");
 972   TString name;
 973
 974   name = "eventAccounting";
 975
 976   // bin index : what it is...
 977   //         0 :  number of event submitted
 978   //         1 :  number accepted by centrality cut
 979   //         2 :  number accepted by centrality cut and z cut
 980   //         3 :  total number of particles that satisfied filter 1
 981   //         4 :  total number of particles that satisfied filter 2
 982   _eventAccounting      = createHisto1D(name,name,10, -0.5, 9.5, "event Code", _title_counts);
 983
 984   name = "m0"; _m0      = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m0, _title_counts);
 985   name = "m1"; _m1      = createHisto1D(name,name,_nBins_M1, _min_M1, _max_M1, _title_m1, _title_counts);
 986   name = "m2"; _m2      = createHisto1D(name,name,_nBins_M2, _min_M2, _max_M2, _title_m2, _title_counts);
 987   name = "m3"; _m3      = createHisto1D(name,name,_nBins_M3, _min_M3, _max_M3, _title_m3, _title_counts);
 988   name = "m4"; _m4      = createHisto1D(name,name,_nBins_M4, _min_M4, _max_M4, _title_m4, _title_counts);
 989   name = "m5"; _m5      = createHisto1D(name,name,_nBins_M5, _min_M5, _max_M5, _title_m5, _title_counts);
 990   name = "m6"; _m6      = createHisto1D(name,name,_nBins_M6, _min_M6, _max_M6, _title_m6, _title_counts);
 991   name = "zV"; _vertexZ = createHisto1D(name,name,_nBins_vertexZ, _min_vertexZ, _max_vertexZ, "z-Vertex (cm)", _title_counts);
 992
 993
 994   name = "Eta";     _etadis   = createHisto1F(name,name, 200, -1.0, 1.0, "#eta","counts");
 995   name = "Phi";     _phidis   = createHisto1F(name,name, 360, 0.0, 6.4, "#phi","counts");
 996   name = "DCAz";    _dcaz     = createHisto1F(name,name, 500, -5.0, 5.0, "dcaZ","counts");
 997   name = "DCAxy";   _dcaxy    = createHisto1F(name,name, 500, -5.0, 5.0, "dcaXY","counts");
 998
 999   /*name = "dedxVsP_1"; _dedxVsP_1  = createHisto2D(name,name,1000,-10.,10.,1000,0.,1000.,"p (GeV/c)", "dedx", "counts");
1000   name = "dedxVsP_2"; _dedxVsP_2  = createHisto2D(name,name,1000,-10.,10.,1000,0.,1000.,"p (GeV/c)", "dedx", "counts");
1001   name = "corrDedxVsP_1"; _corrDedxVsP_1 = createHisto2D(name,name,1000,-10.,10.,1000,0.,500,"p (GeV/c)", "dedx", "counts");
1002   name = "corrDedxVsP_2"; _corrDedxVsP_2 = createHisto2D(name,name,1000,-10.,10.,1000,0.,500,"p (GeV/c)", "dedx", "counts");
1003   */
1004   name = "Nclus1";   _Ncluster1    = createHisto1F(name,name, 200, 0, 200, "Ncluster1","counts");
1005   name = "Nclus2";   _Ncluster2    = createHisto1F(name,name, 200, 0, 200, "Ncluster2","counts");
1006
1007   if (_singlesOnly)
1008     {
1009     name = n1Name+part_1_Name+vsPt;              _n1_1_vsPt              = createHisto1F(name,name, _nBins_pt_1,  _min_pt_1,  _max_pt_1,   _title_pt_1,  _title_AvgN_1);
1010     name = n1Name+part_1_Name+vsZ+vsEtaPhi+vsPt; _n1_1_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_etaPhi_1,  _title_pt_1);
1011     //name = "dedxVsP_1";                          _dedxVsP_1              = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1012     //name = "corrDedxVsP_1";                      _corrDedxVsP_1          = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1013     //name = "betaVsP_1";                          _betaVsP_1              = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
1014
1015     name = n1Name+part_2_Name+vsPt;              _n1_2_vsPt              = createHisto1F(name,name, _nBins_pt_2,  _min_pt_2,  _max_pt_2,   _title_pt_2,  _title_AvgN_2);
1016     name = n1Name+part_2_Name+vsZ+vsEtaPhi+vsPt; _n1_2_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_etaPhi_2,  _title_pt_2);
1017     //name = "dedxVsP_2";                          _dedxVsP_2              = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1018     //name = "corrDedxVsP_2";                      _corrDedxVsP_2          = createHisto2F(name,name,400,-2.,2.,120,0.,120.,"p (GeV/c)", "dedx", "counts");
1019     //name = "betaVsP_2";                          _betaVsP_2              = createHisto2F(name,name,400,-2.,2.,120,0.5,1.1,"p (GeV/c)", "beta", "counts");
1020
1021     }
1022   else
1023     {
1024     name = n1Name+part_1_Name+vsEtaPhi;       _n1_1_vsEtaVsPhi      = createHisto2F(name,name, _nBins_eta_1, _min_eta_1, _max_eta_1,  _nBins_phi_1, _min_phi_1, _max_phi_1,  _title_eta_1,  _title_phi_1,  _title_AvgN_1);
1025     name = s1ptName+part_1_Name+vsEtaPhi;     _s1pt_1_vsEtaVsPhi    = createHisto2F(name,name, _nBins_eta_1, _min_eta_1, _max_eta_1,  _nBins_phi_1, _min_phi_1, _max_phi_1,  _title_eta_1,  _title_phi_1,  _title_AvgSumPt_1);
1026     name = n1Name+part_1_Name+vsM;            _n1_1_vsM             = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
1027     name = s1ptName+part_1_Name+vsM;          _s1pt_1_vsM           = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
1028     name = n1NwName+part_1_Name+vsM;          _n1Nw_1_vsM           = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
1029     name = s1ptNwName+part_1_Name+vsM;        _s1ptNw_1_vsM         = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
1030
1031     name = n1Name+part_2_Name+vsEtaPhi;       _n1_2_vsEtaVsPhi      = createHisto2F(name,name, _nBins_eta_2, _min_eta_2, _max_eta_2,  _nBins_phi_2, _min_phi_2, _max_phi_2,  _title_eta_2,  _title_phi_2,  _title_AvgN_2);
1032     name = s1ptName+part_2_Name+vsEtaPhi;     _s1pt_2_vsEtaVsPhi    = createHisto2F(name,name, _nBins_eta_2, _min_eta_2, _max_eta_2,  _nBins_phi_2, _min_phi_2, _max_phi_2,  _title_eta_2,  _title_phi_2,  _title_AvgSumPt_2);
1033     name = n1Name+part_2_Name + vsM;          _n1_2_vsM             = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_2);
1034     name = s1ptName+part_2_Name + vsM;        _s1pt_2_vsM           = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_2);
1035     name = n1NwName+part_2_Name+vsM;          _n1Nw_2_vsM           = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN_1);
1036     name = s1ptNwName+part_2_Name+vsM;        _s1ptNw_2_vsM         = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPt_1);
1037
1038     name = n2Name+pair_12_Name+vsEtaPhi;      _n2_12_vsEtaPhi       = createHisto1F(name,name, _nBins_etaPhi_12, 0.,        double(_nBins_etaPhi_12), _title_etaPhi_12, _title_AvgN2_12);
1039     name = s2PtPtName+pair_12_Name + vsEtaPhi;_s2PtPt_12_vsEtaPhi   = createHisto1F(name,name, _nBins_etaPhi_12, 0.,        double(_nBins_etaPhi_12), _title_etaPhi_12,  _title_AvgSumPtPt_12);
1040     name = s2PtNName+pair_12_Name + vsEtaPhi; _s2PtN_12_vsEtaPhi    = createHisto1F(name,name, _nBins_etaPhi_12, 0.,        double(_nBins_etaPhi_12), _title_etaPhi_12,  _title_AvgSumPtN_12);
1041     name = s2NPtName+pair_12_Name + vsEtaPhi; _s2NPt_12_vsEtaPhi    = createHisto1F(name,name, _nBins_etaPhi_12, 0.,        double(_nBins_etaPhi_12), _title_etaPhi_12,  _title_AvgNSumPt_12);
1042     name = n2Name+pair_12_Name+vsPtVsPt;      _n2_12_vsPtVsPt       = createHisto2F(name,name, _nBins_pt_1, _min_pt_1, _max_pt_1, _nBins_pt_2, _min_pt_2, _max_pt_2, _title_pt_1, _title_pt_2, _title_AvgN2_12);
1043
1044     name = n2Name+pair_12_Name + vsM;         _n2_12_vsM            = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN2_12);
1045     name = s2PtPtName+pair_12_Name + vsM;     _s2PtPt_12_vsM        = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtPt_12);
1046     name = s2PtNName+pair_12_Name + vsM;      _s2PtN_12_vsM         = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
1047     name = s2NPtName+pair_12_Name + vsM;      _s2NPt_12_vsM         = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
1048
1049     name = n2NwName+pair_12_Name + vsM;       _n2Nw_12_vsM          = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgN2_12);
1050     name = s2PtPtNwName+pair_12_Name + vsM;   _s2PtPtNw_12_vsM      = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtPt_12);
1051     name = s2PtNNwName+pair_12_Name + vsM;    _s2PtNNw_12_vsM       = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgSumPtN_12);
1052     name = s2NPtNwName+pair_12_Name + vsM;    _s2NPtNw_12_vsM       = createProfile(name,name, _nBins_M4, _min_M4, _max_M4, _title_m4, _title_AvgNSumPt_12);
1053
1054     name = "mInv";     _invMass     = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
1055     name = "mInvElec"; _invMassElec = createHisto1F(name,name, 500, 0., 1.000, "M_{inv}","counts");
1056     }
1057
1058   AliInfo(" AliAnalysisTaskDptDptCorrelations::createHistoHistograms() All Done");
1059 }
1060 //-----------------------//
1061
1062 void  AliAnalysisTaskDptDptCorrelations::finalizeHistograms()
1063 {
1064
1065   AliInfo("AliAnalysisTaskDptDptCorrelations::finalizeHistograms() starting");
1066   AliInfo(Form("CorrelationAnalyzers::finalizeHistograms()   _eventCount : %d",int(_eventCount)));
1067   if (_singlesOnly)
1068     {
1069     if (_sameFilter)
1070       {
1071       fillHistoWithArray(_n1_1_vsPt,              __n1_1_vsPt,        _nBins_pt_1);
1072       fillHistoWithArray(_n1_1_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
1073       fillHistoWithArray(_n1_2_vsPt,              __n1_1_vsPt,        _nBins_pt_1);
1074       fillHistoWithArray(_n1_2_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
1075       }
1076     else
1077       {
1078       fillHistoWithArray(_n1_1_vsPt,              __n1_1_vsPt,        _nBins_pt_1);
1079       fillHistoWithArray(_n1_1_vsZVsEtaVsPhiVsPt, __n1_1_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_1, _nBins_pt_1);
1080       fillHistoWithArray(_n1_2_vsPt,              __n1_2_vsPt,        _nBins_pt_2);
1081       fillHistoWithArray(_n1_2_vsZVsEtaVsPhiVsPt, __n1_2_vsZEtaPhiPt, _nBins_vertexZ, _nBins_etaPhi_2, _nBins_pt_2);
1082       }
1083     }
1084   else
1085     {
1086     if (_sameFilter)
1087       {
1088       fillHistoWithArray(_n1_1_vsEtaVsPhi,        __n1_1_vsEtaPhi,    _nBins_eta_1,   _nBins_phi_1);
1089       fillHistoWithArray(_s1pt_1_vsEtaVsPhi,      __s1pt_1_vsEtaPhi,  _nBins_eta_1,   _nBins_phi_1);
1090       fillHistoWithArray(_n1_2_vsEtaVsPhi,        __n1_1_vsEtaPhi,    _nBins_eta_1,   _nBins_phi_1);
1091       fillHistoWithArray(_s1pt_2_vsEtaVsPhi,      __s1pt_1_vsEtaPhi,  _nBins_eta_1,   _nBins_phi_1);
1092       }
1093     else
1094       {
1095       fillHistoWithArray(_n1_1_vsEtaVsPhi,        __n1_1_vsEtaPhi,    _nBins_eta_1,   _nBins_phi_1);
1096       fillHistoWithArray(_s1pt_1_vsEtaVsPhi,      __s1pt_1_vsEtaPhi,  _nBins_eta_1,   _nBins_phi_1);
1097       fillHistoWithArray(_n1_2_vsEtaVsPhi,        __n1_2_vsEtaPhi,    _nBins_eta_2,   _nBins_phi_2);
1098       fillHistoWithArray(_s1pt_2_vsEtaVsPhi,      __s1pt_2_vsEtaPhi,  _nBins_eta_2,   _nBins_phi_2);
1099       }
1100     fillHistoWithArray(_n2_12_vsEtaPhi,     __n2_12_vsEtaPhi,     _nBins_etaPhi_12);
1101     fillHistoWithArray(_s2PtPt_12_vsEtaPhi, __s2ptpt_12_vsEtaPhi, _nBins_etaPhi_12);
1102     fillHistoWithArray(_s2PtN_12_vsEtaPhi,  __s2PtN_12_vsEtaPhi,  _nBins_etaPhi_12);
1103     fillHistoWithArray(_s2NPt_12_vsEtaPhi,  __s2NPt_12_vsEtaPhi,  _nBins_etaPhi_12);
1104     fillHistoWithArray(_n2_12_vsPtVsPt,     __n2_12_vsPtPt,       _nBins_pt_1,    _nBins_pt_2);
1105
1106     }
1107   AliInfo("AliAnalysisTaskDptDptCorrelations::finalizeHistograms()  Done ");
1108 }
1109 //--------------//
1110
1111
1112 void  AliAnalysisTaskDptDptCorrelations::UserExec(Option_t */*option*/)
1113 {
1114
1115   int    k1,k2;
1116   int    iPhi, iEta, iEtaPhi, iPt, charge;
1117   float  q, phi, pt, eta, corr, corrPt, px, py, pz;
1118   int    ij;
1119   //float p, dedx, dedx_1, dedx_2;
1120   //double nSigma;
1121   int    id_1, q_1, iEtaPhi_1, iPt_1;
1122   float  pt_1, px_1, py_1, pz_1, corr_1;
1123   int    id_2, q_2, iEtaPhi_2, iPt_2;
1124   float  pt_2, px_2, py_2, pz_2, corr_2;
1125   float  ptpt;
1126   int    iVertex, iVertexP1, iVertexP2;
1127   int    iZEtaPhiPt;
1128   float  massElecSq = 2.5e-7;
1129   const  AliAODVertex*  vertex;
1130   int    nClus;
1131   bool   bitOK;
1132
1133   AliAnalysisManager* manager = AliAnalysisManager::GetAnalysisManager();
1134   if (!manager) {
1135     //cout<<"ERROR: Analysis manager not found."<<endl;
1136     return;
1137   }
1138   //coneect to the inputHandler------------
1139   AliAODInputHandler* inputHandler = dynamic_cast<AliAODInputHandler*> (manager->GetInputEventHandler());
1140   if (!inputHandler) {
1141     //cout<<"ERROR: Input handler not found."<<endl;
1142     return;
1143   }
1144
1145   fAODEvent = dynamic_cast<AliAODEvent*>(InputEvent());
1146   if (!fAODEvent)
1147     {
1148       //cout<< "ERROR 01: AOD not found " <<endl;
1149       return;
1150     }
1151   fPIDResponse =inputHandler->GetPIDResponse();
1152   if (!fPIDResponse){
1153     AliFatal("This Task needs the PID response attached to the inputHandler");
1154     return;
1155   }
1156
1157   // count all events looked at here
1158   _eventCount++;
1159
1160   // A. Adare - Fix from above to avoid coverity complaint about null dereference
1161   if (_eventAccounting)
1162     {
1163       _eventAccounting->Fill(0);// count all calls to this function
1164     }
1165   else
1166     {
1167
1168     return;
1169     }
1170
1171
1172   //fPIDResponse = ((AliInputEventHandler*)(AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler()))->GetPIDResponse();
1173 //fAODEvent = dynamic_cast<AliAODEvent*> (InputEvent());
1174
1175
1176   _eventAccounting->Fill(1);// count all calls to this function with a valid pointer
1177   //reset single particle counters
1178   k1 = k2 = 0;
1179   __n1_1 = __n1_2 = __s1pt_1 = __s1pt_2 = __n1Nw_1 = __n1Nw_2 = __s1ptNw_1 = __s1ptNw_2 = 0;
1180
1181   float v0Centr  = -999.;
1182   float v0ACentr  = -999.;
1183   float trkCentr = -999.;
1184   float spdCentr = -999.;
1185
1186   float vertexX  = -999;
1187   float vertexY  = -999;
1188   float vertexZ  = -999;
1189   float vertexXY = -999;
1190   //float dcaX     = -999;
1191   //float dcaY     = -999;
1192   float dcaZ     = -999;
1193   float dcaXY    = -999;
1194   float centrality = -999;
1195   //Double_t nSigma =-999;
1196
1197   if(fAODEvent)
1198     {
1199     //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 5" << endl;
1200
1201     //Centrality
1202     AliCentrality* centralityObject =  fAODEvent->GetHeader()->GetCentralityP();
1203     if (centralityObject)
1204       {
1205       //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 6" << endl;
1206
1207       v0Centr  = centralityObject->GetCentralityPercentile("V0M");
1208       v0ACentr  = centralityObject->GetCentralityPercentile("V0A");
1209       trkCentr = centralityObject->GetCentralityPercentile("TRK");
1210       spdCentr = centralityObject->GetCentralityPercentile("CL1");
1211       //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 7" << endl;
1212
1213       }
1214     //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 8" << endl;
1215
1216     // _nTracks  = fAODEvent->GetNTracks(); //OLD
1217     _nTracks  =fAODEvent->GetNumberOfTracks();//NEW Test
1218
1219
1220
1221     //_mult0    = 0;
1222     //_mult1    = 0;
1223     //_mult2    = 0;
1224     _mult3    = _nTracks;
1225     _mult4    = v0Centr;
1226     _mult4a    = v0ACentr;
1227     _mult5    = trkCentr;
1228     _mult6    = spdCentr;
1229     _field    = fAODEvent->GetMagneticField();
1230
1231     //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - _field:" << _field << endl;
1232
1233     //_centralityMethod
1234     switch (_centralityMethod)
1235       {
1236       case 0: centrality = _mult0; break;
1237       case 1: centrality = _mult1; break;
1238       case 2: centrality = _mult2; break;
1239       case 3: centrality = _mult3; break;
1240       case 4: centrality = _mult4; break;
1241       case 5: centrality = _mult5; break;
1242       case 6: centrality = _mult6; break;
1243       case 7: centrality = _mult4a; break;
1244       }
1245
1246
1247     if ( centrality < _centralityMin ||
1248          centrality > _centralityMax ||
1249          fabs(v0Centr-trkCentr)>5.0)
1250       {
1251
1252       return;
1253       }
1254     _eventAccounting->Fill(2);// count all events with right centrality
1255     //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 12" << endl;
1256
1257     // filter on z and xy vertex
1258     vertex = (AliAODVertex*) fAODEvent->GetPrimaryVertexSPD();
1259     if (!vertex || vertex->GetNContributors()<1)
1260       {
1261       vertexZ  = -999;
1262       vertexXY = -999;
1263       AliInfo("AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - No valid vertex object or poor vertex");
1264       }
1265     else
1266       {
1267         vertexX = vertex->GetX();
1268         vertexY = vertex->GetY();
1269         vertexZ = vertex->GetZ();
1270         vertexXY = sqrt(vertexX*vertexX+vertexY*vertexY);
1271       }
1272     if (!vertex ||
1273         vertexZ  < _vertexZMin  ||
1274         vertexZ  > _vertexZMax  ||
1275         vertexXY < _vertexXYMin ||
1276         vertexXY > _vertexXYMax)
1277       return;
1278
1279     iVertex = int((vertexZ-_min_vertexZ)/_width_vertexZ);
1280     iVertexP1 = iVertex*_nBins_etaPhiPt_1;
1281     iVertexP2 = iVertex*_nBins_etaPhiPt_2;
1282     if (iVertex<0 || iVertex>=_nBins_vertexZ)
1283       {
1284       AliError("AliAnalysisTaskDptDptCorrelations::Exec(Option_t *option) iVertex<0 || iVertex>=_nBins_vertexZ ");
1285       return;
1286       }
1287     _eventAccounting->Fill(3);// count all calls to this function with a valid pointer
1288     //======================
1289
1290     TExMap *trackMap = new TExMap();//Mapping matrix----
1291
1292     //1st loop track
1293     for(Int_t i = 0; i < _nTracks; i++)
1294       {
1295         AliAODTrack* aodTrack = dynamic_cast<AliAODTrack *>(fAODEvent->GetTrack(i));
1296         if(!aodTrack) {
1297           AliError(Form("ERROR: Could not retrieve AODtrack %d",i));
1298           continue;
1299         }
1300         //Double_t tpcSignalAll = aodTrack->GetTPCsignal();
1301         //_dedxVsP_1->Fill(aodTrack->GetTPCmomentum(),tpcSignalAll);
1302         Int_t gID = aodTrack->GetID();
1303         if (aodTrack->TestFilterBit(1)) trackMap->Add(gID, i);//Global tracks
1304       }
1305
1306
1307     //--------------------------
1308
1309     AliAODTrack* newAodTrack;
1310
1311     //2nd loop track starts here
1312     for (int iTrack=0; iTrack< _nTracks; iTrack++)
1313       {
1314       //cout << "AliAnalysisTaskDptDptCorrelations::UserExec(Option_t *option) - 16 iTrack:" << iTrack << endl;
1315
1316       AliAODTrack * t = (AliAODTrack *) fAODEvent->GetTrack(iTrack);
1317       if (!t)
1318         {
1319         AliError(Form("AliAnalysisTaskDptDptCorrelations::Exec(Option_t *option) No track ofr iTrack=%d", iTrack));
1320         continue;
1321         }
1322
1323       bitOK  = t->TestFilterBit(_trackFilterBit);
1324       if (!bitOK) continue;
1325
1326       Int_t gID = t->GetID();
1327       newAodTrack = gID >= 0 ?t : fAODEvent->GetTrack(trackMap->GetValue(-1-gID));
1328
1329
1330       q      = t->Charge();
1331       charge = int(q);
1332       phi    = t->Phi();
1333       //p      = t->P();
1334       pt     = t->Pt();
1335       px     = t->Px();
1336       py     = t->Py();
1337       pz     = t->Pz();
1338       eta    = t->Eta();
1339       //dedx   = t->GetTPCsignal();
1340       nClus  = t->GetTPCNcls();
1341       Double_t nclus2 = t->GetTPCClusterInfo(2,1);
1342
1343
1344
1345
1346       //Float_t nsigmaTPCPID = -999.;
1347       //Float_t nsigmaTOFPID = -999.;
1348
1349       Double_t nsigmaelectron = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kElectron));
1350       Double_t nsigmapion = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kPion));
1351       Double_t nsigmakaon = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kKaon));
1352       Double_t nsigmaproton = TMath::Abs(fPIDResponse->NumberOfSigmasTPC(newAodTrack,(AliPID::EParticleType)AliPID::kProton));
1353
1354       //Fill QA before the PID
1355       //_dedxVsP_1 -> Fill(newAodTrack->P()*newAodTrack->Charge(),newAodTrack->GetTPCsignal());
1356       //_corrDedxVsP_1 -> Fill(newAodTrack->P()*newAodTrack->Charge(),nsigmaelectron);
1357       //end of QA-before pid
1358
1359       //-----------------------------------------
1360       //nsigma cut to reject electron (NEW)
1361       if(nsigmaelectron  < fNSigmaCut
1362          && nsigmapion   > fNSigmaCut
1363          && nsigmakaon   > fNSigmaCut
1364          && nsigmaproton > fNSigmaCut ) continue;
1365       //------------------------------------------
1366
1367       //Fill QA after the PID
1368       // _dedxVsP_2 -> Fill(newAodTrack->P()*newAodTrack->Charge(),newAodTrack->GetTPCsignal());
1369       //_corrDedxVsP_2 -> Fill(newAodTrack->P()*newAodTrack->Charge(),nsigmaelectron);
1370       //-----------------
1371
1372       //Particle 1
1373       if (_requestedCharge_1 == charge &&
1374           pt       >=  _min_pt_1 &&
1375           pt       <   _max_pt_1 &&
1376           eta      >=  _min_eta_1 &&
1377           eta      <   _max_eta_1)
1378         {
1379
1380
1381           iPhi   = int( phi/_width_phi_1);
1382
1383           if (iPhi<0 || iPhi>=_nBins_phi_1 )
1384             {
1385               AliWarning("AliAnalysisTaskDptDptCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1386           return;
1387           }
1388
1389         iEta    = int((eta-_min_eta_1)/_width_eta_1);
1390         if (iEta<0 || iEta>=_nBins_eta_1)
1391           {
1392           AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1393           continue;
1394           }
1395         iPt     = int((pt -_min_pt_1 )/_width_pt_1 );
1396         if (iPt<0  || iPt >=_nBins_pt_1)
1397           {
1398           AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1399           continue;
1400           }
1401         iEtaPhi = iEta*_nBins_phi_1+iPhi;
1402         iZEtaPhiPt = iVertexP1 + iEtaPhi*_nBins_pt_1 + iPt;
1403
1404         if (_correctionWeight_1)
1405           corr = _correctionWeight_1[iZEtaPhiPt];
1406         else
1407           corr = 1;
1408         if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1)
1409           {
1410           AliWarning("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_1");
1411           continue;
1412           }
1413
1414
1415         if (_singlesOnly)
1416           {
1417
1418           __n1_1_vsPt[iPt]               += corr;          //cout << "step 15" << endl;
1419           __n1_1_vsZEtaPhiPt[iZEtaPhiPt] += corr;       //cout << "step 12" << endl;
1420
1421           }
1422         else
1423           {
1424           corrPt                      = corr*pt;
1425           _id_1[k1]                   = iTrack;
1426           _charge_1[k1]               = charge;
1427           _iEtaPhi_1[k1]              = iEtaPhi;
1428           _iPt_1[k1]                  = iPt;
1429           _pt_1[k1]                   = pt;
1430           _px_1[k1]                   = px;
1431           _py_1[k1]                   = py;
1432           _pz_1[k1]                   = pz;
1433           _correction_1[k1]           = corr;
1434           __n1_1                      += corr;
1435           __n1_1_vsEtaPhi[iEtaPhi]    += corr;
1436           __s1pt_1                    += corrPt;
1437           __s1pt_1_vsEtaPhi[iEtaPhi]  += corrPt;
1438           __n1Nw_1                    += 1;
1439           __s1ptNw_1                  += pt;
1440           ++k1;
1441           if (k1>=arraySize)
1442             {
1443             AliError(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) k1 >=arraySize; arraySize: %d",arraySize));
1444             return;
1445             }
1446           }
1447         }
1448
1449       if (!_sameFilter &&
1450           _requestedCharge_2 == charge &&
1451           pt       >=  _min_pt_2 &&
1452           pt       <   _max_pt_2 &&
1453           eta      >=  _min_eta_2 &&
1454           eta      <   _max_eta_2)
1455         {
1456
1457           dcaXY = t->DCA(); //new change Prabhat
1458           dcaZ  = t->ZAtDCA(); //new change Prabhat
1459
1460           //Check dca cuts for systematics
1461           //Comment this portion when not required
1462           /*if (dcaZ  <  _dcaZMin ||
1463               dcaZ  >  _dcaZMax ||
1464               dcaXY <  _dcaXYMin ||
1465               dcaXY >  _dcaXYMax ) continue;
1466           */
1467           //==== QA ===========================
1468           _dcaz->Fill(dcaZ);
1469           _dcaxy->Fill(dcaXY);
1470           _etadis->Fill(eta);
1471           _phidis->Fill(phi);
1472           _Ncluster1->Fill(nClus);
1473           _Ncluster2->Fill(nclus2);
1474
1475           //===================================
1476
1477           iPhi   = int( phi/_width_phi_2);
1478
1479         if (iPhi<0 || iPhi>=_nBins_phi_2 )
1480           {
1481           AliWarning("AliAnalysisTaskDptDptCorrelations::analyze() iPhi<0 || iPhi>=_nBins_phi_1");
1482           return;
1483           }
1484
1485         iEta    = int((eta-_min_eta_2)/_width_eta_2);
1486         if (iEta<0 || iEta>=_nBins_eta_2)
1487           {
1488           AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iEta: %d", iEta));
1489           continue;
1490           }
1491         iPt     = int((pt -_min_pt_2 )/_width_pt_2 );
1492         if (iPt<0  || iPt >=_nBins_pt_2)
1493           {
1494           AliWarning(Form("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) Mismatched iPt: %d",iPt));
1495           continue;
1496           }
1497
1498         iEtaPhi = iEta*_nBins_phi_2+iPhi;
1499         iZEtaPhiPt = iVertexP2 + iEtaPhi*_nBins_pt_2 + iPt;
1500         if (iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2)
1501           {
1502           AliWarning("AliAnalysisTaskDptDptCorrelations::analyze(AliceEvent * event) iZEtaPhiPt<0 || iZEtaPhiPt>=_nBins_zEtaPhiPt_2");
1503           continue;
1504           }
1505
1506
1507         if (_correctionWeight_2)
1508           corr = _correctionWeight_2[iZEtaPhiPt];
1509         else
1510           corr = 1;
1511         //dpt = pt - (charge>0) ? _avgPt_vsEtaPhi_2p[iEtaPhi] : _avgPt_vsEtaPhi_2m[iEtaPhi];
1512
1513         if (_singlesOnly)
1514           {
1515           //_dedxVsP_2->Fill(p*q,dedx);
1516           __n1_2_vsPt[iPt]               += corr;          //cout << "step 15" << endl;
1517           __n1_2_vsZEtaPhiPt[iZEtaPhiPt] += corr;       //cout << "step 12" << endl;
1518           }
1519         else
1520           {
1521           corrPt                      = corr*pt;
1522           _id_2[k2]                   = iTrack;         //cout << "step 1" << endl;
1523           _charge_2[k2]               = charge;         //cout << "step 2" << endl;
1524           _iEtaPhi_2[k2]              = iEtaPhi;        //cout << "step 3" << endl;
1525           _iPt_2[k2]                  = iPt;            //cout << "step 4" << endl;
1526           _pt_2[k2]                   = pt;             //cout << "step 5" << endl;
1527           _px_2[k2]                   = px;             //cout << "step 6" << endl;
1528           _py_2[k2]                   = py;             //cout << "step 7" << endl;
1529           _pz_2[k2]                   = pz;             //cout << "step 8" << endl;
1530           _correction_2[k2]           = corr;           //cout << "step 9" << endl;
1531           __n1_2                      += corr;          //cout << "step 10" << endl;
1532           __s1pt_2                    += corrPt;        //cout << "step 13" << endl;
1533           __n1Nw_2                    += 1;
1534           __n1_2_vsEtaPhi[iEtaPhi]    += corr;          //cout << "step 11" << endl;
1535           __s1pt_2_vsEtaPhi[iEtaPhi]  += corrPt;        //cout << "step 14" << endl;
1536           __s1ptNw_2                  += pt;
1537           ++k2;
1538           if (k2>=arraySize)
1539             {
1540             AliWarning(Form("-W- k2 >=arraySize; arraySize: %d",arraySize));
1541             return;
1542             }
1543           }
1544
1545       //cout << "done with track" << endl;
1546       } //iTrack
1547     } //aod
1548     }
1549
1550
1551   //cout << "Filling histograms now" << endl;
1552   _m0->Fill(_mult0);
1553   _m1->Fill(_mult1);
1554   _m2->Fill(_mult2);
1555   _m3->Fill(_mult3);
1556   _m4->Fill(_mult4);
1557   _m5->Fill(_mult5);
1558   _m6->Fill(_mult6);
1559   _vertexZ->Fill(vertexZ);
1560
1561   if (_singlesOnly)
1562     {
1563     // nothing to do here.
1564     }
1565   else
1566     {
1567     if (_sameFilter)
1568       {
1569       _n1_1_vsM->Fill(centrality,      __n1_1);
1570       _s1pt_1_vsM->Fill(centrality,    __s1pt_1);
1571       _n1Nw_1_vsM->Fill(centrality,    __n1Nw_1);
1572       _s1ptNw_1_vsM->Fill(centrality,  __s1ptNw_1);
1573       _n1_2_vsM->Fill(centrality,      __n1_1);
1574       _s1pt_2_vsM->Fill(centrality,    __s1pt_1);
1575       _n1Nw_2_vsM->Fill(centrality,    __n1Nw_1);
1576       _s1ptNw_2_vsM->Fill(centrality,  __s1ptNw_1);
1577       // reset pair counters
1578       __n2_12   = __s2ptpt_12   = __s2NPt_12    = __s2PtN_12    = 0;
1579       __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12  = __s2PtNNw_12  = 0;
1580       if (_field>0)
1581         {
1582         for (int i1=0; i1<k1; i1++)
1583           {
1584           ////cout << "         i1:" << i1 << endl;
1585           id_1      = _id_1[i1];           ////cout << "       id_1:" << id_1 << endl;
1586           q_1       = _charge_1[i1];       ////cout << "        q_1:" << q_1 << endl;
1587           iEtaPhi_1 = _iEtaPhi_1[i1];      ////cout << "  iEtaPhi_1:" << iEtaPhi_1 << endl;
1588           iPt_1     = _iPt_1[i1];          ////cout << "      iPt_1:" << iPt_1 << endl;
1589           corr_1    = _correction_1[i1];   ////cout << "     corr_1:" << corr_1 << endl;
1590           pt_1      = _pt_1[i1];           ////cout << "       pt_1:" << pt_1 << endl;
1591           //dedx_1    = _dedx_1[i1];         ////cout << "     dedx_1:" << dedx_1 << endl;
1592                                            //1 and 2
1593           for (int i2=i1+1; i2<k1; i2++)
1594             {
1595               ////cout << "         i2:" << i2 << endl;
1596               id_2      = _id_1[i2];              ////cout << "       id_2:" << id_2 << endl;
1597               if (id_1!=id_2)
1598                 {
1599                 q_2       = _charge_1[i2];     ////cout << "        q_1:" << q_1 << endl;
1600                 iEtaPhi_2 = _iEtaPhi_1[i2];    ////cout << "  iEtaPhi_1:" << iEtaPhi_1 << endl;
1601                 iPt_2     = _iPt_1[i2];        ////cout << "      iPt_1:" << iPt_1 << endl;
1602                 corr_2    = _correction_1[i2]; ////cout << "     corr_1:" << corr_1 << endl;
1603                 pt_2      = _pt_1[i2];         ////cout << "       pt_1:" << pt_1 << endl;
1604                 //dedx_2    = _dedx_1[i2];       ////cout << "     dedx_2:" << dedx_2 << endl;
1605                 corr      = corr_1*corr_2;
1606                 if (q_2>q_1 || (q_1>0 && q_2>0 && pt_2<=pt_1) || (q_1<0 && q_2<0 && pt_2>=pt_1))
1607                   {
1608                   ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2;   ////cout << " ij:" << ij<< endl;
1609                   }
1610                 else // swap particles
1611                   {
1612                   ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1;   ////cout << " ij:" << ij<< endl;
1613                   }
1614
1615                 __n2_12                  += corr;
1616                 __n2_12_vsEtaPhi[ij]     += corr;
1617                 ptpt                     = pt_1*pt_2;
1618                 __s2ptpt_12              += corr*ptpt;
1619                 __s2PtN_12               += corr*pt_1;
1620                 __s2NPt_12               += corr*pt_2;
1621                 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1622                 __s2PtN_12_vsEtaPhi[ij]  += corr*pt_1;
1623                 __s2NPt_12_vsEtaPhi[ij]  += corr*pt_2;
1624                 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1625
1626                 __n2Nw_12                  += 1;
1627                 __s2ptptNw_12              += ptpt;
1628                 __s2PtNNw_12               += pt_1;
1629                 __s2NPtNw_12               += pt_2;
1630
1631                 }
1632             } //i2
1633           } //i1
1634         }
1635       else // field<0
1636         {
1637         for (int i1=0; i1<k1; i1++)
1638           {
1639           ////cout << "         i1:" << i1 << endl;
1640           id_1      = _id_1[i1];           ////cout << "       id_1:" << id_1 << endl;
1641           q_1       = _charge_1[i1];       ////cout << "        q_1:" << q_1 << endl;
1642           iEtaPhi_1 = _iEtaPhi_1[i1];      ////cout << "  iEtaPhi_1:" << iEtaPhi_1 << endl;
1643           iPt_1     = _iPt_1[i1];          ////cout << "      iPt_1:" << iPt_1 << endl;
1644           corr_1    = _correction_1[i1];   ////cout << "     corr_1:" << corr_1 << endl;
1645           pt_1      = _pt_1[i1];           ////cout << "       pt_1:" << pt_1 << endl;
1646           //dedx_1    = _dedx_1[i1];         ////cout << "     dedx_1:" << dedx_1 << endl;
1647                                            //1 and 2
1648           for (int i2=i1+1; i2<k1; i2++)
1649             {
1650               ////cout << "         i2:" << i2 << endl;
1651               id_2      = _id_1[i2];              ////cout << "       id_2:" << id_2 << endl;
1652               if (id_1!=id_2)
1653                 {
1654                 q_2       = _charge_1[i2];     ////cout << "        q_2:" << q_2 << endl;
1655                 iEtaPhi_2 = _iEtaPhi_1[i2];    ////cout << "  iEtaPhi_2:" << iEtaPhi_2 << endl;
1656                 iPt_2     = _iPt_1[i2];        ////cout << "      iPt_2:" << iPt_2 << endl;
1657                 corr_2    = _correction_1[i2]; ////cout << "     corr_2:" << corr_2 << endl;
1658                 pt_2      = _pt_1[i2];         ////cout << "       pt_2:" << pt_2 << endl;
1659                 //dedx_2    = _dedx_1[i2];       ////cout << "     dedx_2:" << dedx_2 << endl;
1660                 corr      = corr_1*corr_2;
1661                 if ( q_2<q_1 || (q_1>0 && q_2>0 && pt_2>=pt_1) || (q_1<0 && q_2<0 && pt_2<=pt_1))
1662                   {
1663                   ij = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2;   ////cout << " ij:" << ij<< endl;
1664                   }
1665                 else // swap particles
1666                   {
1667                   ij = iEtaPhi_2*_nBins_etaPhi_1 + iEtaPhi_1;   ////cout << " ij:" << ij<< endl;
1668                   }
1669
1670                 __n2_12                  += corr;
1671                 __n2_12_vsEtaPhi[ij]     += corr;
1672                 ptpt                     = pt_1*pt_2;
1673                 __s2ptpt_12              += corr*ptpt;
1674                 __s2PtN_12               += corr*pt_1;
1675                 __s2NPt_12               += corr*pt_2;
1676                 __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1677                 __s2PtN_12_vsEtaPhi[ij]  += corr*pt_1;
1678                 __s2NPt_12_vsEtaPhi[ij]  += corr*pt_2;
1679                 __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1680
1681                 __n2Nw_12                  += 1;
1682                 __s2ptptNw_12              += ptpt;
1683                 __s2PtNNw_12               += pt_1;
1684                 __s2NPtNw_12               += pt_2;
1685
1686                 }
1687             } //i2
1688           } //i1
1689         }
1690       }
1691     else  // filter 1 and 2 are different -- must do all particle pairs...
1692       {
1693       _n1_1_vsM->Fill(centrality,      __n1_1);
1694       _s1pt_1_vsM->Fill(centrality,    __s1pt_1);
1695       _n1Nw_1_vsM->Fill(centrality,    __n1Nw_1);
1696       _s1ptNw_1_vsM->Fill(centrality,  __s1ptNw_1);
1697       _n1_2_vsM->Fill(centrality,      __n1_2);
1698       _s1pt_2_vsM->Fill(centrality,    __s1pt_2);
1699       _n1Nw_2_vsM->Fill(centrality,    __n1Nw_2);
1700       _s1ptNw_2_vsM->Fill(centrality,  __s1ptNw_2);
1701       // reset pair counters
1702       __n2_12   = __s2ptpt_12   = __s2NPt_12    = __s2PtN_12    = 0;
1703       __n2Nw_12 = __s2ptptNw_12 = __s2NPtNw_12  = __s2PtNNw_12  = 0;
1704       for (int i1=0; i1<k1; i1++)
1705         {
1706         ////cout << "         i1:" << i1 << endl;
1707         id_1      = _id_1[i1];           ////cout << "       id_1:" << id_1 << endl;
1708         q_1       = _charge_1[i1];       ////cout << "        q_1:" << q_1 << endl;
1709         iEtaPhi_1 = _iEtaPhi_1[i1];      ////cout << "  iEtaPhi_1:" << iEtaPhi_1 << endl;
1710         iPt_1     = _iPt_1[i1];          ////cout << "      iPt_1:" << iPt_1 << endl;
1711         corr_1    = _correction_1[i1];   ////cout << "     corr_1:" << corr_1 << endl;
1712         pt_1      = _pt_1[i1];           ////cout << "       pt_1:" << pt_1 << endl;
1713         px_1      = _px_1[i1];          ////cout << "      px_1:" << px_1 << endl;
1714         py_1      = _py_1[i1];          ////cout << "      py_1:" << py_1 << endl;
1715         pz_1      = _pz_1[i1];          ////cout << "      pz_1:" << pz_1 << endl;
1716         //dedx_1    = _dedx_1[i1];        ////cout << "     dedx_1:" << dedx_1 << endl;
1717
1718         //1 and 2
1719         for (int i2=0; i2<k2; i2++)
1720           {
1721             ////cout << "         i2:" << i2 << endl;
1722             id_2   = _id_2[i2];              ////cout << "       id_2:" << id_2 << endl;
1723             if (id_1!=id_2)  // exclude auto correlation
1724               {
1725               q_2       = _charge_2[i2];     ////cout << "        q_2:" << q_2 << endl;
1726               iEtaPhi_2 = _iEtaPhi_2[i2];    ////cout << "  iEtaPhi_2:" << iEtaPhi_2 << endl;
1727               iPt_2     = _iPt_2[i2];        ////cout << "      iPt_2:" << iPt_2 << endl;
1728               corr_2    = _correction_2[i2]; ////cout << "     corr_2:" << corr_2 << endl;
1729               pt_2      = _pt_2[i2];         ////cout << "       pt_2:" << pt_2 << endl;
1730               px_2      = _px_2[i2];          ////cout << "      px_2:" << px_2 << endl;
1731               py_2      = _py_2[i2];          ////cout << "      py_2:" << py_2 << endl;
1732               pz_2      = _pz_2[i2];          ////cout << "      pz_2:" << pz_2 << endl;
1733               //dedx_2    = _dedx_2[i2];        ////cout << "     dedx_2:" << dedx_2 << endl;
1734
1735
1736               if (_rejectPairConversion)
1737                 {
1738                 float e1Sq = massElecSq + pt_1*pt_1 + pz_1*pz_1;
1739                 float e2Sq = massElecSq + pt_2*pt_2 + pz_2*pz_2;
1740                 float mInvSq = 2*(massElecSq + sqrt(e1Sq*e2Sq) - px_1*px_2 - py_1*py_2 - pz_1*pz_2 );
1741                 float mInv = sqrt(mInvSq);
1742                 _invMass->Fill(mInv);
1743                 //if (mInv<0.5)
1744                 //{
1745                     // if (dedx_1>75. && dedx_2>75.)
1746                     //{
1747                     _invMassElec->Fill(mInv);
1748                     //if (mInv<0.05) continue;
1749                     //}
1750                     //}
1751                 }
1752               corr      = corr_1*corr_2;
1753               ij        = iEtaPhi_1*_nBins_etaPhi_1 + iEtaPhi_2;   ////cout << " ij:" << ij<< endl;
1754               __n2_12                  += corr;
1755               __n2_12_vsEtaPhi[ij]     += corr;
1756               ptpt                     = pt_1*pt_2;
1757               __s2ptpt_12              += corr*ptpt;
1758               __s2PtN_12               += corr*pt_1;
1759               __s2NPt_12               += corr*pt_2;
1760               __s2ptpt_12_vsEtaPhi[ij] += corr*ptpt;
1761               __s2PtN_12_vsEtaPhi[ij]  += corr*pt_1;
1762               __s2NPt_12_vsEtaPhi[ij]  += corr*pt_2;
1763               __n2_12_vsPtPt[iPt_1*_nBins_pt_2 + iPt_2] += corr;
1764               __n2Nw_12                  += 1;
1765               __s2ptptNw_12              += ptpt;
1766               __s2PtNNw_12               += pt_1;
1767               __s2NPtNw_12               += pt_2;
1768
1769               }
1770           } //i2
1771         } //i1
1772       }
1773
1774     _n2_12_vsM->Fill(centrality,     __n2_12);
1775     _s2PtPt_12_vsM->Fill(centrality, __s2ptpt_12);
1776     _s2PtN_12_vsM->Fill(centrality,  __s2NPt_12);
1777     _s2NPt_12_vsM->Fill(centrality,  __s2PtN_12);
1778
1779     _n2Nw_12_vsM->Fill(centrality,     __n2Nw_12);
1780     _s2PtPtNw_12_vsM->Fill(centrality, __s2ptptNw_12);
1781     _s2PtNNw_12_vsM->Fill(centrality,  __s2NPtNw_12);
1782     _s2NPtNw_12_vsM->Fill(centrality,  __s2PtNNw_12);
1783
1784     }
1785
1786
1787   AliInfo("AliAnalysisTaskDptDptCorrelations::UserExec()   -----------------Event Done ");
1788   PostData(0,_outputHistoList);
1789
1790 }
1791
1792 void   AliAnalysisTaskDptDptCorrelations::FinishTaskOutput()
1793 {
1794   AliInfo("AliAnalysisTaskDptDptCorrelations::FinishTaskOutput() Starting.");
1795   Printf("= 0 ====================================================================");
1796   finalizeHistograms();
1797   AliInfo("= 1 ====================================================================");
1798   PostData(0,_outputHistoList);
1799   AliInfo("= 2 ====================================================================");
1800   AliInfo("AliAnalysisTaskDptDptCorrelations::FinishTaskOutput() Done.");
1801 }
1802
1803 void   AliAnalysisTaskDptDptCorrelations::Terminate(Option_t* /*option*/)
1804 {
1805   AliInfo("AliAnalysisTaskDptDptCorrelations::Terminate() Starting/Done.");
1806 }
1807
1808
1809 //Tools
1810 //===================================================================================================
1811 void  AliAnalysisTaskDptDptCorrelations::fillHistoWithArray(TH1 * h, float * array, int size)
1812 {
1813   int i, i1;
1814   float v1, ev1, v2, ev2, sum, esum;
1815   for (i=0, i1=1; i<size; ++i,++i1)
1816     {
1817     v1  = array[i]; ev1 = sqrt(v1);
1818     v2  = h->GetBinContent(i1);
1819     ev2 = h->GetBinError(i1);
1820     sum = v1 + v2;
1821     esum = sqrt(ev1*ev1+ev2*ev2);
1822     h->SetBinContent(i1,sum);
1823     h->SetBinError(i1,esum);
1824     }
1825 }
1826
1827 void  AliAnalysisTaskDptDptCorrelations::fillHistoWithArray(TH2 * h, float * array, int size1, int size2)
1828 {
1829   int i, i1;
1830   int j, j1;
1831   float v1, ev1, v2, ev2, sum, esum;
1832   for (i=0, i1=1; i<size1; ++i,++i1)
1833     {
1834     for (j=0, j1=1; j<size2; ++j,++j1)
1835       {
1836       v1  = array[i*size2+j]; ev1 = sqrt(v1);
1837       v2  = h->GetBinContent(i1,j1);
1838       ev2 = h->GetBinError(i1,j1);
1839       sum = v1 + v2;
1840       esum = sqrt(ev1*ev1+ev2*ev2);
1841       h->SetBinContent(i1,j1,sum);
1842       h->SetBinError(i1,j1,esum);
1843       }
1844     }
1845 }
1846
1847 void  AliAnalysisTaskDptDptCorrelations::fillHistoWithArray(TH3 * h, float * array, int size1, int size2, int size3)
1848 {
1849   int i, i1;
1850   int j, j1;
1851   int k, k1;
1852   float v1, ev1, v2, ev2, sum, esum;
1853   int size23 = size2*size3;
1854   for (i=0, i1=1; i<size1; ++i,++i1)
1855     {
1856     for (j=0, j1=1; j<size2; ++j,++j1)
1857       {
1858       for (k=0, k1=1; k<size3; ++k,++k1)
1859         {
1860         v1  = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
1861         v2  = h->GetBinContent(i1,j1,k1);
1862         ev2 = h->GetBinError(i1,j1,k1);
1863         sum = v1 + v2;
1864         esum = sqrt(ev1*ev1+ev2*ev2);
1865         h->SetBinContent(i1,j1,k1,sum);
1866         h->SetBinError(i1,j1,k1,esum);
1867         }
1868       }
1869     }
1870 }
1871
1872 void  AliAnalysisTaskDptDptCorrelations::fillHistoWithArray(TH1 * h, double * array, int size)
1873 {
1874   int i, i1;
1875   double v1, ev1, v2, ev2, sum, esum;
1876   for (i=0, i1=1; i<size; ++i,++i1)
1877     {
1878     v1  = array[i]; ev1 = sqrt(v1);
1879     v2  = h->GetBinContent(i1);
1880     ev2 = h->GetBinError(i1);
1881     sum = v1 + v2;
1882     esum = sqrt(ev1*ev1+ev2*ev2);
1883     h->SetBinContent(i1,sum);
1884     h->SetBinError(i1,esum);
1885     }
1886 }
1887
1888 void  AliAnalysisTaskDptDptCorrelations::fillHistoWithArray(TH2 * h, double * array, int size1, int size2)
1889 {
1890   int i, i1;
1891   int j, j1;
1892   double v1, ev1, v2, ev2, sum, esum;
1893   for (i=0, i1=1; i<size1; ++i,++i1)
1894     {
1895     for (j=0, j1=1; j<size2; ++j,++j1)
1896       {
1897       v1  = array[i*size2+j]; ev1 = sqrt(v1);
1898       v2  = h->GetBinContent(i1,j1);
1899       ev2 = h->GetBinError(i1,j1);
1900       sum = v1 + v2;
1901       esum = sqrt(ev1*ev1+ev2*ev2);
1902       h->SetBinContent(i1,j1,sum);
1903       h->SetBinError(i1,j1,esum);
1904       }
1905     }
1906 }
1907
1908 void  AliAnalysisTaskDptDptCorrelations::fillHistoWithArray(TH3 * h, double * array, int size1, int size2, int size3)
1909 {
1910   int i, i1;
1911   int j, j1;
1912   int k, k1;
1913   double v1, ev1, v2, ev2, sum, esum;
1914   int size23 = size2*size3;
1915   for (i=0, i1=1; i<size1; ++i,++i1)
1916     {
1917     for (j=0, j1=1; j<size2; ++j,++j1)
1918       {
1919       for (k=0, k1=1; k<size3; ++k,++k1)
1920         {
1921         v1  = array[i*size23+j*size3+k]; ev1 = sqrt(v1);
1922         v2  = h->GetBinContent(i1,j1,k1);
1923         ev2 = h->GetBinError(i1,j1,k1);
1924         sum = v1 + v2;
1925         esum = sqrt(ev1*ev1+ev2*ev2);
1926         h->SetBinContent(i1,j1,k1,sum);
1927         h->SetBinError(i1,j1,k1,esum);
1928         }
1929       }
1930     }
1931 }
1932
1933 //________________________________________________________________________
1934 double *  AliAnalysisTaskDptDptCorrelations::getDoubleArray(int size, double v)
1935 {
1936   /// Allocate an array of type double with n values
1937   /// Initialize the array to the given value
1938   double * array = new double [size];
1939   for (int i=0;i<size;++i) array[i]=v;
1940   return array;
1941 }
1942
1943 //________________________________________________________________________
1944 float *  AliAnalysisTaskDptDptCorrelations::getFloatArray(int size, float v)
1945 {
1946   /// Allocate an array of type float with n values
1947   /// Initialize the array to the given value
1948   float * array = new float [size];
1949   for (int i=0;i<size;++i) array[i]=v;
1950   return array;
1951 }
1952
1953
1954 //________________________________________________________________________
1955 TH1D * AliAnalysisTaskDptDptCorrelations::createHisto1D(const TString &  name, const TString &  title,
1956                                                       int n, double xMin, double xMax,
1957                                                       const TString &  xTitle, const TString &  yTitle)
1958 {
1959   //CreateHisto new 1D historgram
1960   AliInfo(Form("createHisto 1D histo %s   nBins: %d  xMin: %f   xMax: %f",name.Data(),n,xMin,xMax));
1961   TH1D * h = new TH1D(name,title,n,xMin,xMax);
1962   h->GetXaxis()->SetTitle(xTitle);
1963   h->GetYaxis()->SetTitle(yTitle);
1964   addToList(h);
1965   return h;
1966 }
1967
1968
1969 //________________________________________________________________________
1970 TH1D * AliAnalysisTaskDptDptCorrelations::createHisto1D(const TString &  name, const TString &  title,
1971                                                       int n, double * bins,
1972                                                       const TString &  xTitle, const TString &  yTitle)
1973 {
1974   AliInfo(Form("createHisto 1D histo %s   with %d non uniform nBins",name.Data(),n));
1975   TH1D * h = new TH1D(name,title,n,bins);
1976   h->GetXaxis()->SetTitle(xTitle);
1977   h->GetYaxis()->SetTitle(yTitle);
1978   addToList(h);
1979   return h;
1980 }
1981
1982
1983 //________________________________________________________________________
1984 TH2D * AliAnalysisTaskDptDptCorrelations::createHisto2D(const TString &  name, const TString &  title,
1985                                                       int nx, double xMin, double xMax, int ny, double yMin, double yMax,
1986                                                       const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
1987 {
1988   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));
1989   TH2D * h = new TH2D(name,title,nx,xMin,xMax,ny,yMin,yMax);
1990   h->GetXaxis()->SetTitle(xTitle);
1991   h->GetYaxis()->SetTitle(yTitle);
1992   h->GetZaxis()->SetTitle(zTitle);
1993   addToList(h);
1994   return h;
1995 }
1996
1997 //________________________________________________________________________
1998 TH2D * AliAnalysisTaskDptDptCorrelations::createHisto2D(const TString &  name, const TString &  title,
1999                                                       int nx, double* xbins, int ny, double yMin, double yMax,
2000                                                       const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
2001 {
2002   AliInfo(Form("createHisto 2D histo %s   with %d non uniform nBins",name.Data(),nx));
2003   TH2D * h;
2004   h = new TH2D(name,title,nx,xbins,ny,yMin,yMax);
2005   h->GetXaxis()->SetTitle(xTitle);
2006   h->GetYaxis()->SetTitle(yTitle);
2007   h->GetZaxis()->SetTitle(zTitle);
2008   addToList(h);
2009   return h;
2010 }
2011
2012 //// F /////
2013 //________________________________________________________________________
2014 TH1F * AliAnalysisTaskDptDptCorrelations::createHisto1F(const TString &  name, const TString &  title,
2015                                                         int n, double xMin, double xMax,
2016                                                         const TString &  xTitle, const TString &  yTitle)
2017 {
2018   //CreateHisto new 1D historgram
2019   AliInfo(Form("createHisto 1D histo %s   nBins: %d  xMin: %f   xMax: %f",name.Data(),n,xMin,xMax));
2020   TH1F * h = new TH1F(name,title,n,xMin,xMax);
2021   h->GetXaxis()->SetTitle(xTitle);
2022   h->GetYaxis()->SetTitle(yTitle);
2023   addToList(h);
2024   return h;
2025 }
2026
2027
2028 //________________________________________________________________________
2029 TH1F * AliAnalysisTaskDptDptCorrelations::createHisto1F(const TString &  name, const TString &  title,
2030                                                         int n, double * bins,
2031                                                         const TString &  xTitle, const TString &  yTitle)
2032 {
2033   AliInfo(Form("createHisto 1D histo %s   with %d non uniform nBins",name.Data(),n));
2034   TH1F * h = new TH1F(name,title,n,bins);
2035   h->GetXaxis()->SetTitle(xTitle);
2036   h->GetYaxis()->SetTitle(yTitle);
2037   addToList(h);
2038   return h;
2039 }
2040
2041
2042 //________________________________________________________________________
2043 TH2F * AliAnalysisTaskDptDptCorrelations::createHisto2F(const TString &  name, const TString &  title,
2044                                                         int nx, double xMin, double xMax, int ny, double yMin, double yMax,
2045                                                         const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
2046 {
2047   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));
2048   TH2F * h = new TH2F(name,title,nx,xMin,xMax,ny,yMin,yMax);
2049   h->GetXaxis()->SetTitle(xTitle);
2050   h->GetYaxis()->SetTitle(yTitle);
2051   h->GetZaxis()->SetTitle(zTitle);
2052   addToList(h);
2053   return h;
2054 }
2055
2056 //________________________________________________________________________
2057 TH2F * AliAnalysisTaskDptDptCorrelations::createHisto2F(const TString &  name, const TString &  title,
2058                                                         int nx, double* xbins, int ny, double yMin, double yMax,
2059                                                         const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
2060 {
2061   AliInfo(Form("createHisto 2D histo %s   with %d non uniform nBins",name.Data(),nx));
2062   TH2F * h;
2063   h = new TH2F(name,title,nx,xbins,ny,yMin,yMax);
2064   h->GetXaxis()->SetTitle(xTitle);
2065   h->GetYaxis()->SetTitle(yTitle);
2066   h->GetZaxis()->SetTitle(zTitle);
2067   addToList(h);
2068   return h;
2069 }
2070
2071
2072 //________________________________________________________________________
2073 TH3F * AliAnalysisTaskDptDptCorrelations::createHisto3F(const TString &  name, const TString &  title,
2074                                                       int nx, double xMin, double xMax,
2075                                                       int ny, double yMin, double yMax,
2076                                                       int nz, double zMin, double zMax,
2077                                                       const TString &  xTitle, const TString &  yTitle, const TString &  zTitle)
2078 {
2079   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));
2080   TH3F * h = new TH3F(name,title,nx,xMin,xMax,ny,yMin,yMax,nz,zMin,zMax);
2081   h->GetXaxis()->SetTitle(xTitle);
2082   h->GetYaxis()->SetTitle(yTitle);
2083   h->GetZaxis()->SetTitle(zTitle);
2084   addToList(h);
2085   return h;
2086 }
2087
2088
2089 //________________________________________________________________________
2090 TProfile * AliAnalysisTaskDptDptCorrelations::createProfile(const TString & name, const TString & description,
2091                                                             int nx,double xMin,double xMax,
2092                                                             const TString &  xTitle, const TString &  yTitle)
2093 {
2094   AliInfo(Form("createHisto 1D profile %s   nBins: %d  xMin: %f10.4 xMax: %f10.4",name.Data(),nx,xMin,xMax));
2095   TProfile * h = new TProfile(name,description,nx,xMin,xMax);
2096   h->GetXaxis()->SetTitle(xTitle);
2097   h->GetYaxis()->SetTitle(yTitle);
2098   addToList(h);
2099   return h;
2100 }
2101
2102 //________________________________________________________________________
2103 TProfile * AliAnalysisTaskDptDptCorrelations::createProfile(const TString &  name,const TString &  description,
2104                                                             int nx,  double* bins,
2105                                                             const TString &  xTitle, const TString &  yTitle)
2106 {
2107   AliInfo(Form("createHisto 1D profile %s  with %d non uniform bins",name.Data(),nx));
2108   TProfile * h = new TProfile(name,description,nx,bins);
2109   h->GetXaxis()->SetTitle(xTitle);
2110   h->GetYaxis()->SetTitle(yTitle);
2111   addToList(h);
2112   return h;
2113 }
2114
2115
2116 void   AliAnalysisTaskDptDptCorrelations::addToList(TH1 *h)
2117 {
2118   if (_outputHistoList)
2119     {
2120     _outputHistoList->Add(h);
2121     }
2122   else
2123     AliInfo("addToList(TH1 *h) _outputHistoList is null!!!!! Should abort ship");
2124
2125 }
2126
2127
2128