PWGCF/Correlations/DPhi/AliAnalysisTaskpzpz.cxx

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