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