PWGLF/SPECTRA/ChargedHadrons/dNdPt/AlidNdPtAnalysispPb.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15 //------------------------------------------------------------------------------
  16 // AlidNdPtAnalysispPb class.
  17 //
  18 // a. functionality:
  19 // - fills analysis control histograms
  20 // - fills generic correction matrices
  21 // - generates correction matrices
  22 //
  23 // b. data members:
  24 // - generic correction matrices
  25 // - control histograms
  26 //
  27 // last change: 2013-06-19 by M.Knichel
  28 //
  29 // meaning of different multiplicities:
  30 // multRec      : number of reconstructed tracks, after AcceptanceCuts and TrackCuts
  31 // multRecMult  : number of reconstructed tracks, after MultAcceptanceCuts and MultTrackCuts
  32 // multMB       : number of contributers to vertex
  33 // multTrueMC   : MC true mult, after MultAcceptanceCuts
  34 // mutlAll      : number of ESD tracks
  35 // mutlAcc      : number of ESD tracks after AcceptanceCuts
  36 //------------------------------------------------------------------------------
  37
  38 #include "TH1.h"
  39 #include "TH2.h"
  40 #include "TCanvas.h"
  41 #include "THnSparse.h"
  42
  43 #include "AliHeader.h"
  44 #include "AliInputEventHandler.h"
  45 #include "AliAnalysisManager.h"
  46 #include "AliGenEventHeader.h"
  47 #include "AliStack.h"
  48 #include "AliESDEvent.h"
  49 #include "AliMCEvent.h"
  50 #include "AliESDtrackCuts.h"
  51 #include "AliLog.h"
  52 #include "AliMultiplicity.h"
  53 #include "AliTracker.h"
  54
  55 #include "AlidNdPtEventCuts.h"
  56 #include "AlidNdPtAcceptanceCuts.h"
  57 #include "AliPhysicsSelection.h"
  58 #include "AliTriggerAnalysis.h"
  59 #include "AliCentrality.h"
  60 #include "AliAnalysisUtils.h"
  61
  62 #include "AliPWG0Helper.h"
  63 #include "AlidNdPtHelper.h"
  64 #include "AlidNdPtAnalysispPb.h"
  65
  66 using namespace std;
  67
  68 ClassImp(AlidNdPtAnalysispPb)
  69
  70 //_____________________________________________________________________________
  71   AlidNdPtAnalysispPb::AlidNdPtAnalysispPb(): AlidNdPt(),
  72   fAnalysisFolder(0),
  73   fHistogramsOn(kFALSE),
  74
  75   // event multiplicity correlation matrix
  76   fEventMultCorrelationMatrix(0),
  77
  78   // rec. track pt vs true track pt correlation matrix
  79   fTrackPtCorrelationMatrix(0),
  80
  81   // event level correction
  82   fGenEventMatrix(0),
  83   fGenSDEventMatrix(0),
  84   fGenDDEventMatrix(0),
  85   fGenNDEventMatrix(0),
  86   fGenNSDEventMatrix(0),
  87
  88   fTriggerEventMatrix(0),
  89   fTriggerSDEventMatrix(0),
  90   fTriggerDDEventMatrix(0),
  91   fTriggerNDEventMatrix(0),
  92   fTriggerNSDEventMatrix(0),
  93
  94   fRecEventMatrix(0),
  95   fRecSDEventMatrix(0),
  96   fRecDDEventMatrix(0),
  97   fRecNDEventMatrix(0),
  98   fRecNSDEventMatrix(0),
  99
 100   //
 101   // track-event level correction
 102   //
 103   fGenTrackEventMatrix(0),
 104   fGenTrackSDEventMatrix(0),
 105   fGenTrackDDEventMatrix(0),
 106   fGenTrackNDEventMatrix(0),
 107   fGenTrackNSDEventMatrix(0),
 108
 109   fTriggerTrackEventMatrix(0),
 110   fTriggerTrackSDEventMatrix(0),
 111   fTriggerTrackDDEventMatrix(0),
 112   fTriggerTrackNDEventMatrix(0),
 113   fTriggerTrackNSDEventMatrix(0),
 114
 115   fRecTrackEventMatrix(0),
 116   fRecTrackSDEventMatrix(0),
 117   fRecTrackDDEventMatrix(0),
 118   fRecTrackNDEventMatrix(0),
 119   fRecTrackNSDEventMatrix(0),
 120
 121   // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
 122   fGenTrackMatrix(0),
 123   fGenPrimTrackMatrix(0),
 124   fRecPrimTrackMatrix(0),
 125
 126   // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)
 127   fRecTrackMatrix(0),
 128   fRecSecTrackMatrix(0),
 129
 130   // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix)
 131   fRecMultTrackMatrix(0),
 132
 133   // event control histograms
 134   fMCEventHist1(0),
 135   fRecEventHist1(0),
 136   fRecEventHist2(0),
 137   fRecMCEventHist1(0),
 138   fRecMCEventHist2(0),
 139   fRecMCEventHist3(0),
 140
 141   // rec. pt and eta resolution w.r.t MC
 142   fRecMCTrackHist1(0),
 143
 144   //multple reconstructed tracks
 145   fMCMultRecTrackHist1(0),
 146
 147   // rec. track control histograms
 148   fRecTrackHist2(0),
 149
 150   // Generic histograms to be corrected
 151   fRecEventHist(0),
 152   fRecTrackHist(0),
 153   fEventCount(0),
 154   fEventMultHist(0),
 155   fMCPrimTrackHist(0),
 156
 157   // Candle event histogram
 158   fRecCandleEventMatrix(0),
 159
 160   fCentralityEventHist(0),
 161   fCentralityTrackHist(0),
 162   fCentralityEstimatorsList(0),
 163   fDimensionsCentralityEstimators(0),
 164   fCentralityNbins(0),
 165   fCentralityNedges(0),
 166   fBinsCentrality(0),
 167   fNVCentralityEvent(0),
 168   fNVCentralityTrack(0),
 169   fVCentralityEvent(0),
 170   fVCentralityTrack(0),
 171
 172   fMultNbins(0),
 173   fMultNbinsTE(0),
 174   fPtNbins(0),
 175   fPtCorrNbins(0),
 176   fEtaNbins(0),
 177   fZvNbins(0),
 178   fMultNedges(0),
 179   fMultNedgesTE(0),
 180   fPtNedges(0),
 181   fPtCorrNedges(0),
 182   fEtaNedges(0),
 183   fZvNedges(0),
 184   fBinsMult(0),
 185   fBinsMultTE(0),
 186   fBinsPt(0),
 187   fBinsPtCorr(0),
 188   fBinsEta(0),
 189   fBinsZv(0),
 190
 191   fRapidityShift(-4.65409416218532379e-01),
 192   fUtils(0),
 193   fIs2013pA(kTRUE),
 194
 195   fIsInit(kFALSE)
 196
 197 {
 198   // default constructor
 199   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
 200     fMCTrackHist1[i]=0;
 201     fMCPrimTrackHist1[i]=0;
 202     fMCPrimTrackHist2[i]=0;
 203     fMCSecTrackHist1[i]=0;
 204     fRecTrackHist1[i]=0;
 205     fRecTrackMultHist1[i]=0;
 206   }
 207   //Init();
 208 }
 209
 210 //_____________________________________________________________________________
 211 AlidNdPtAnalysispPb::AlidNdPtAnalysispPb(Char_t* name, Char_t* title): AlidNdPt(name,title),
 212   fAnalysisFolder(0),
 213   fHistogramsOn(kFALSE),
 214
 215   // event multiplicity correlation matrix
 216   fEventMultCorrelationMatrix(0),
 217
 218   // rec. track pt vs true track pt correlation matrix
 219   fTrackPtCorrelationMatrix(0),
 220
 221   // event level correction
 222   fGenEventMatrix(0),
 223   fGenSDEventMatrix(0),
 224   fGenDDEventMatrix(0),
 225   fGenNDEventMatrix(0),
 226   fGenNSDEventMatrix(0),
 227
 228   fTriggerEventMatrix(0),
 229   fTriggerSDEventMatrix(0),
 230   fTriggerDDEventMatrix(0),
 231   fTriggerNDEventMatrix(0),
 232   fTriggerNSDEventMatrix(0),
 233
 234   fRecEventMatrix(0),
 235   fRecSDEventMatrix(0),
 236   fRecDDEventMatrix(0),
 237   fRecNDEventMatrix(0),
 238   fRecNSDEventMatrix(0),
 239
 240   //
 241   // track-event level correction
 242   //
 243   fGenTrackEventMatrix(0),
 244   fGenTrackSDEventMatrix(0),
 245   fGenTrackDDEventMatrix(0),
 246   fGenTrackNDEventMatrix(0),
 247   fGenTrackNSDEventMatrix(0),
 248
 249   fTriggerTrackEventMatrix(0),
 250   fTriggerTrackSDEventMatrix(0),
 251   fTriggerTrackDDEventMatrix(0),
 252   fTriggerTrackNDEventMatrix(0),
 253   fTriggerTrackNSDEventMatrix(0),
 254
 255   fRecTrackEventMatrix(0),
 256   fRecTrackSDEventMatrix(0),
 257   fRecTrackDDEventMatrix(0),
 258   fRecTrackNDEventMatrix(0),
 259   fRecTrackNSDEventMatrix(0),
 260
 261   // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
 262   fGenTrackMatrix(0),
 263   fGenPrimTrackMatrix(0),
 264   fRecPrimTrackMatrix(0),
 265
 266   // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)
 267   fRecTrackMatrix(0),
 268   fRecSecTrackMatrix(0),
 269
 270   // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix)
 271   fRecMultTrackMatrix(0),
 272
 273   // event control histograms
 274   fMCEventHist1(0),
 275   fRecEventHist1(0),
 276   fRecEventHist2(0),
 277   fRecMCEventHist1(0),
 278   fRecMCEventHist2(0),
 279   fRecMCEventHist3(0),
 280
 281   // rec. pt and eta resolution w.r.t MC
 282   fRecMCTrackHist1(0),
 283
 284   //multple reconstructed tracks
 285   fMCMultRecTrackHist1(0),
 286
 287   // rec. track control histograms
 288   fRecTrackHist2(0),
 289
 290   // Generic histograms to be corrected
 291   fRecEventHist(0),
 292   fRecTrackHist(0),
 293   fEventCount(0),
 294   fEventMultHist(0),
 295   fMCPrimTrackHist(0),
 296
 297   // Candle event histogram
 298   fRecCandleEventMatrix(0),
 299
 300   fCentralityEventHist(0),
 301   fCentralityTrackHist(0),
 302   fCentralityEstimatorsList(0),
 303   fDimensionsCentralityEstimators(0),
 304   fCentralityNbins(0),
 305   fCentralityNedges(0),
 306   fBinsCentrality(0),
 307   fNVCentralityEvent(0),
 308   fNVCentralityTrack(0),
 309   fVCentralityEvent(0),
 310   fVCentralityTrack(0),
 311
 312   fMultNbins(0),
 313   fMultNbinsTE(0),
 314   fPtNbins(0),
 315   fPtCorrNbins(0),
 316   fEtaNbins(0),
 317   fZvNbins(0),
 318   fMultNedges(0),
 319   fMultNedgesTE(0),
 320   fPtNedges(0),
 321   fPtCorrNedges(0),
 322   fEtaNedges(0),
 323   fZvNedges(0),
 324   fBinsMult(0),
 325   fBinsMultTE(0),
 326   fBinsPt(0),
 327   fBinsPtCorr(0),
 328   fBinsEta(0),
 329   fBinsZv(0),
 330
 331   fRapidityShift(-4.65409416218532379e-01),
 332   fUtils(0),
 333   fIs2013pA(kTRUE),
 334
 335   fIsInit(kFALSE)
 336
 337 {
 338   //
 339   // constructor
 340   //
 341   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
 342     fMCTrackHist1[i]=0;
 343     fMCPrimTrackHist1[i]=0;
 344     fMCPrimTrackHist2[i]=0;
 345     fMCSecTrackHist1[i]=0;
 346     fRecTrackHist1[i]=0;
 347     fRecTrackMultHist1[i]=0;
 348   }
 349
 350   //Init();
 351 }
 352
 353 //_____________________________________________________________________________
 354 AlidNdPtAnalysispPb::~AlidNdPtAnalysispPb() {
 355   //
 356   // destructor
 357   //
 358   if(fEventMultCorrelationMatrix) delete fEventMultCorrelationMatrix; fEventMultCorrelationMatrix=0;
 359   //
 360   if(fTrackPtCorrelationMatrix) delete fTrackPtCorrelationMatrix; fTrackPtCorrelationMatrix=0;
 361   //
 362   if(fGenEventMatrix) delete fGenEventMatrix; fGenEventMatrix=0;
 363   if(fGenSDEventMatrix) delete fGenSDEventMatrix; fGenSDEventMatrix=0;
 364   if(fGenDDEventMatrix) delete fGenDDEventMatrix; fGenDDEventMatrix=0;
 365   if(fGenNDEventMatrix) delete fGenNDEventMatrix; fGenNDEventMatrix=0;
 366   if(fGenNSDEventMatrix) delete fGenNSDEventMatrix; fGenNSDEventMatrix=0;
 367
 368   if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;
 369   if(fTriggerSDEventMatrix) delete fTriggerSDEventMatrix; fTriggerSDEventMatrix=0;
 370   if(fTriggerDDEventMatrix) delete fTriggerDDEventMatrix; fTriggerDDEventMatrix=0;
 371   if(fTriggerNDEventMatrix) delete fTriggerNDEventMatrix; fTriggerNDEventMatrix=0;
 372   if(fTriggerNSDEventMatrix) delete fTriggerNSDEventMatrix; fTriggerNSDEventMatrix=0;
 373
 374   if(fRecEventMatrix) delete fRecEventMatrix; fRecEventMatrix=0;
 375   if(fRecSDEventMatrix) delete fRecSDEventMatrix; fRecSDEventMatrix=0;
 376   if(fRecDDEventMatrix) delete fRecDDEventMatrix; fRecDDEventMatrix=0;
 377   if(fRecNDEventMatrix) delete fRecNDEventMatrix; fRecNDEventMatrix=0;
 378   if(fRecNSDEventMatrix) delete fRecNSDEventMatrix; fRecNSDEventMatrix=0;
 379
 380   if(fRecCandleEventMatrix) delete fRecCandleEventMatrix; fRecCandleEventMatrix=0;
 381   //
 382   if(fGenTrackEventMatrix) delete fGenTrackEventMatrix; fGenTrackEventMatrix=0;
 383   if(fGenTrackSDEventMatrix) delete fGenTrackSDEventMatrix; fGenTrackSDEventMatrix=0;
 384   if(fGenTrackDDEventMatrix) delete fGenTrackDDEventMatrix; fGenTrackDDEventMatrix=0;
 385   if(fGenTrackNDEventMatrix) delete fGenTrackNDEventMatrix; fGenTrackNDEventMatrix=0;
 386   if(fGenTrackNSDEventMatrix) delete fGenTrackNSDEventMatrix; fGenTrackNSDEventMatrix=0;
 387
 388   if(fTriggerTrackEventMatrix) delete fTriggerTrackEventMatrix; fTriggerTrackEventMatrix=0;
 389   if(fTriggerTrackSDEventMatrix) delete fTriggerTrackSDEventMatrix; fTriggerTrackSDEventMatrix=0;
 390   if(fTriggerTrackDDEventMatrix) delete fTriggerTrackDDEventMatrix; fTriggerTrackDDEventMatrix=0;
 391   if(fTriggerTrackNDEventMatrix) delete fTriggerTrackNDEventMatrix; fTriggerTrackNDEventMatrix=0;
 392   if(fTriggerTrackNSDEventMatrix) delete fTriggerTrackNSDEventMatrix; fTriggerTrackNSDEventMatrix=0;
 393
 394   if(fRecTrackEventMatrix) delete fRecTrackEventMatrix; fRecTrackEventMatrix=0;
 395   if(fRecTrackSDEventMatrix) delete fRecTrackSDEventMatrix; fRecTrackSDEventMatrix=0;
 396   if(fRecTrackDDEventMatrix) delete fRecTrackDDEventMatrix; fRecTrackDDEventMatrix=0;
 397   if(fRecTrackNDEventMatrix) delete fRecTrackNDEventMatrix; fRecTrackNDEventMatrix=0;
 398   if(fRecTrackNSDEventMatrix) delete fRecTrackNSDEventMatrix; fRecTrackNSDEventMatrix=0;
 399
 400   //
 401   if(fGenTrackMatrix) delete fGenTrackMatrix; fGenTrackMatrix=0;
 402   if(fGenPrimTrackMatrix) delete fGenPrimTrackMatrix; fGenPrimTrackMatrix=0;
 403   if(fRecPrimTrackMatrix) delete fRecPrimTrackMatrix; fRecPrimTrackMatrix=0;
 404   //
 405   if(fRecTrackMatrix) delete fRecTrackMatrix; fRecTrackMatrix=0;
 406   if(fRecSecTrackMatrix) delete fRecSecTrackMatrix; fRecSecTrackMatrix=0;
 407   //
 408   if(fRecMultTrackMatrix) delete fRecMultTrackMatrix; fRecMultTrackMatrix=0;
 409   //
 410   // Control histograms
 411   //
 412   if(fMCEventHist1) delete fMCEventHist1; fMCEventHist1=0;
 413   if(fRecEventHist1) delete fRecEventHist1; fRecEventHist1=0;
 414   if(fRecEventHist2) delete fRecEventHist2; fRecEventHist2=0;
 415   if(fRecMCEventHist1) delete fRecMCEventHist1; fRecMCEventHist1=0;
 416   if(fRecMCEventHist2) delete fRecMCEventHist2; fRecMCEventHist2=0;
 417   if(fRecMCEventHist3) delete fRecMCEventHist3; fRecMCEventHist3=0;
 418   //
 419   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
 420     if(fMCTrackHist1[i]) delete fMCTrackHist1[i]; fMCTrackHist1[i]=0;
 421     if(fMCPrimTrackHist1[i]) delete fMCPrimTrackHist1[i]; fMCPrimTrackHist1[i]=0;
 422     if(fMCPrimTrackHist2[i]) delete fMCPrimTrackHist2[i]; fMCPrimTrackHist2[i]=0;
 423     if(fMCSecTrackHist1[i]) delete fMCSecTrackHist1[i]; fMCSecTrackHist1[i]=0;
 424     if(fRecTrackHist1[i]) delete fRecTrackHist1[i]; fRecTrackHist1[i]=0;
 425     if(fRecTrackMultHist1[i]) delete fRecTrackMultHist1[i]; fRecTrackMultHist1[i]=0;
 426   }
 427   if(fRecMCTrackHist1) delete fRecMCTrackHist1; fRecMCTrackHist1=0;
 428   if(fMCMultRecTrackHist1) delete fMCMultRecTrackHist1; fMCMultRecTrackHist1=0;
 429   if(fRecTrackHist2) delete fRecTrackHist2; fRecTrackHist2=0;
 430
 431   //
 432   if(fRecEventHist) delete fRecEventHist; fRecEventHist=0;
 433   if(fRecTrackHist) delete fRecTrackHist; fRecTrackHist=0;
 434   if(fEventCount) delete fEventCount; fEventCount=0;
 435   if(fEventMultHist) delete fEventMultHist; fEventMultHist=0;
 436   if(fMCPrimTrackHist) delete fMCPrimTrackHist; fMCPrimTrackHist=0;
 437
 438   //
 439   if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
 440
 441   if (fBinsMult) delete[] fBinsMult; fBinsMult=0;
 442   if (fBinsPt) delete[] fBinsPt; fBinsPt=0;
 443   if (fBinsPtCorr) delete[] fBinsPtCorr; fBinsPtCorr=0;
 444   if (fBinsEta) delete[] fBinsEta; fBinsEta=0;
 445   if (fBinsZv) delete[] fBinsMult; fBinsZv=0;
 446
 447   if (fCentralityEstimatorsList) delete fCentralityEstimatorsList; fCentralityEstimatorsList=0;
 448   if (fBinsCentrality) delete[] fBinsCentrality; fBinsCentrality=0;
 449   if (fCentralityEventHist) delete fCentralityEventHist; fCentralityEventHist=0;
 450   if (fCentralityTrackHist) delete fCentralityTrackHist; fCentralityTrackHist=0;
 451   if (fVCentralityEvent) delete[] fVCentralityEvent; fVCentralityEvent=0;
 452   if (fVCentralityTrack) delete[] fVCentralityTrack; fVCentralityTrack=0;
 453
 454   if (fUtils) delete fUtils; fUtils=0;
 455
 456 }
 457
 458
 459 //_____________________________________________________________________________
 460 void AlidNdPtAnalysispPb::SetCentralityEstimators(const char* estimators)
 461 {
 462   if (CanChangeBins()) {
 463     if (fCentralityEstimatorsList) delete fCentralityEstimatorsList; fCentralityEstimatorsList=0;
 464     TString estimatorsString(estimators);
 465     fCentralityEstimatorsList = estimatorsString.Tokenize(" ,;:");
 466     fDimensionsCentralityEstimators = fCentralityEstimatorsList->GetEntries();
 467   }
 468 }
 469
 470
 471 //_____________________________________________________________________________
 472 Bool_t AlidNdPtAnalysispPb::CanChangeBins()
 473 {
 474   if (fIsInit) {
 475       AliDebug(AliLog::kError, "Object AlidNdPtAnalysispPb already initialized. Cannot change.");
 476       return kFALSE;
 477   }
 478   return kTRUE;
 479 }
 480
 481
 482 //_____________________________________________________________________________
 483 void AlidNdPtAnalysispPb::Init()
 484 {
 485     //define default binning
 486     Double_t binsMultDefault[15] = {-0.5, 0.5 , 1.5 , 2.5 , 3.5 , 4.5 , 5.5 , 6.5 , 7.5 , 8.5, 9.5, 10.5, 20.5, 50.5, 150.5};
 487     Double_t binsPtDefault[69] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,2.8,3.0,3.2,3.4,3.6,3.8,4.0,4.5,5.0,5.5,6.0,6.5,7.0,8.0,9.0,10.0,11.0,12.0,13.0,14.0,15.0,16.0,18.0,20.0,22.0,24.0,26.0,28.0,30.0,32.0,34.0,36.0,40.0,45.0,50.0};
 488     Double_t binsPtCorrDefault[37] = {0.,0.05,0.1,0.15,0.2,0.25,0.3,0.35,0.4,0.45,0.5,0.55,0.6,0.65,0.7,0.75,0.8,0.85,0.9,0.95,1.0,1.1,1.2,1.3,1.4,1.5,1.6,1.7,1.8,1.9,2.0,2.2,2.4,2.6,3.0,4.0,50.0};
 489     Double_t binsEtaDefault[31] = {-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5};
 490     Double_t binsZvDefault[13] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
 491
 492     Double_t binCentralityDefault[6] = {0.,20.,40.,60.,80.,100.};
 493
 494    // if no binning is set, use the default
 495    if (!fBinsMult)        { SetBinsMult(14,binsMultDefault); }
 496    if (!fBinsMultTE)      { SetBinsMultTE(14,binsMultDefault); }
 497    if (!fBinsPt)          { SetBinsPt(68,binsPtDefault); }
 498    if (!fBinsPtCorr)      { SetBinsPtCorr(36,binsPtCorrDefault); }
 499    if (!fBinsEta)         { SetBinsEta(30,binsEtaDefault); }
 500    if (!fBinsZv)          { SetBinsZv(12,binsZvDefault); }
 501    if (!fBinsCentrality)  { SetBinsCentrality(5,binCentralityDefault); }
 502
 503   // Int_t binsTrackMatrix[3]={fZvNbins,fPtCorrNbins,fEtaNbins};
 504   Int_t binsTrackEventCorrMatrix[4]={fZvNbins,fPtCorrNbins,fEtaNbins,fMultNbinsTE};
 505   Int_t binsTrackEventCorrMatrixCent[4]={fZvNbins,fPtCorrNbins,fEtaNbins,fCentralityNbins};
 506
 507
 508   //
 509   // Generic histograms to be corrected
 510   //
 511   Int_t binsEventHist[2]={fZvNbins,fMultNbins};
 512   //Double_t minEventHist[2]={-fBinsZv[0],fBinsMult[0]};
 513   //Double_t maxEventHist[2]={fBinsZv[fZvNbins],fBinsMult[fMultNbins]};
 514
 515   fRecEventHist = new THnSparseF("fRecEventHist","Zv:multMB",2,binsEventHist); //,minEventHist,maxEventHist);
 516   fRecEventHist->SetBinEdges(0,fBinsZv);
 517   fRecEventHist->SetBinEdges(1,fBinsMult);
 518   fRecEventHist->GetAxis(0)->SetTitle("Zv (cm)");
 519   fRecEventHist->GetAxis(1)->SetTitle("multiplicity MB");
 520   fRecEventHist->Sumw2();
 521
 522   //
 523   Int_t binsTrackHist[4]={fZvNbins,fPtNbins,fEtaNbins,fMultNbins};
 524  // Double_t minTrackHist[4]={-25.,0.,-1.5,-0.5};
 525  // Double_t maxTrackHist[4]={25.,50.,1.5,149.5};
 526
 527   fRecTrackHist = new THnSparseF("fRecTrackHist","Zv:pT:eta:multRecMult",4,binsTrackHist); //,minTrackHist,maxTrackHist);
 528   fRecTrackHist->SetBinEdges(0,fBinsZv);
 529   fRecTrackHist->SetBinEdges(1,fBinsPt);
 530   fRecTrackHist->SetBinEdges(2,fBinsEta);
 531   fRecTrackHist->SetBinEdges(3,fBinsMult);
 532   fRecTrackHist->GetAxis(0)->SetTitle("Zv (cm)");
 533   fRecTrackHist->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
 534   fRecTrackHist->GetAxis(2)->SetTitle("#eta");
 535   fRecTrackHist->GetAxis(3)->SetTitle("multiplicity (multCuts)");
 536   fRecTrackHist->Sumw2();
 537
 538
 539   Int_t binsMCPrimTrackHist[4]={fZvNbins,fPtNbins,fEtaNbins,fMultNbins};
 540  // Double_t minTrackHist[4]={-25.,0.,-1.5,-0.5};
 541  // Double_t maxTrackHist[4]={25.,50.,1.5,149.5};
 542
 543   fMCPrimTrackHist = new THnSparseF("fMCPrimTrackHist","Zv:mcpT:mceta:multTrueMC",4,binsMCPrimTrackHist);
 544   fMCPrimTrackHist->SetBinEdges(0,fBinsZv);
 545   fMCPrimTrackHist->SetBinEdges(1,fBinsPt);
 546   fMCPrimTrackHist->SetBinEdges(2,fBinsEta);
 547   fMCPrimTrackHist->SetBinEdges(3,fBinsMult);
 548   fMCPrimTrackHist->GetAxis(0)->SetTitle("Zv (cm)");
 549   fMCPrimTrackHist->GetAxis(1)->SetTitle("MC p_{T} (GeV/c)");
 550   fMCPrimTrackHist->GetAxis(2)->SetTitle("#eta (MC)");
 551   fMCPrimTrackHist->GetAxis(3)->SetTitle("true multiplicity (MC)");
 552   fMCPrimTrackHist->Sumw2();
 553
 554   //
 555   // rec. vs MC correlation matrices
 556   //
 557   Int_t binsMultTrueEventMatrix[3]={fMultNbins,fMultNbins,fMultNbins};
 558 //   Double_t minMultTrueEventMatrix[3]={-0.5,-0.5,-0.5};
 559 //   Double_t maxMultTrueEventMatrix[3]={149.5,149.5,149.5};
 560   fEventMultCorrelationMatrix = new THnSparseF("fEventMultCorrelationMatrix","multRecMult:multTrueMC:multMB",3,binsMultTrueEventMatrix); //,minMultTrueEventMatrix,maxMultTrueEventMatrix);
 561   fEventMultCorrelationMatrix->SetBinEdges(0,fBinsMult);
 562   fEventMultCorrelationMatrix->SetBinEdges(1,fBinsMult);
 563   fEventMultCorrelationMatrix->SetBinEdges(2,fBinsMult);
 564   fEventMultCorrelationMatrix->GetAxis(0)->SetTitle("multiplicity (multCuts)");
 565   fEventMultCorrelationMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 566   fEventMultCorrelationMatrix->GetAxis(2)->SetTitle("MB multiplicity");
 567   fEventMultCorrelationMatrix->Sumw2();
 568
 569   Int_t binsTrackPtCorrelationMatrix[3]={fPtCorrNbins,fPtCorrNbins,fEtaNbins};
 570   fTrackPtCorrelationMatrix = new THnSparseF("fTrackPtCorrelationMatrix","Pt:mcPt:mcEta",3,binsTrackPtCorrelationMatrix);
 571   fTrackPtCorrelationMatrix->SetBinEdges(0,fBinsPtCorr);
 572   fTrackPtCorrelationMatrix->SetBinEdges(1,fBinsPtCorr);
 573   fTrackPtCorrelationMatrix->SetBinEdges(2,fBinsEta);
 574   fTrackPtCorrelationMatrix->GetAxis(0)->SetTitle("Pt (GeV/c)");
 575   fTrackPtCorrelationMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 576   fTrackPtCorrelationMatrix->GetAxis(2)->SetTitle("mcEta");
 577   fTrackPtCorrelationMatrix->Sumw2();
 578
 579   //
 580   // Efficiency and contamination correction matrices
 581   //
 582   Int_t binsEventMatrix[2]={fZvNbins,fMultNbins};
 583   Int_t binsEventMatrixCent[3]={fZvNbins,fMultNbins,fCentralityNbins};
 584 //   Double_t minEventMatrix[2]={-25.,-0.5};
 585 //   Double_t maxEventMatrix[2]={25.,149.5};
 586
 587   fGenEventMatrix = new THnSparseF("fGenEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 588   fGenEventMatrix->SetBinEdges(0,fBinsZv);
 589   fGenEventMatrix->SetBinEdges(1,fBinsMult);
 590   fGenEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 591   fGenEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 592   fGenEventMatrix->Sumw2();
 593
 594   fGenSDEventMatrix = new THnSparseF("fGenSDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 595   fGenSDEventMatrix->SetBinEdges(0,fBinsZv);
 596   fGenSDEventMatrix->SetBinEdges(1,fBinsMult);
 597   fGenSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 598   fGenSDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 599   fGenSDEventMatrix->Sumw2();
 600
 601   fGenDDEventMatrix = new THnSparseF("fGenDDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 602   fGenDDEventMatrix->SetBinEdges(0,fBinsZv);
 603   fGenDDEventMatrix->SetBinEdges(1,fBinsMult);
 604   fGenDDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 605   fGenDDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 606   fGenDDEventMatrix->Sumw2();
 607
 608   fGenNDEventMatrix = new THnSparseF("fGenNDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 609   fGenNDEventMatrix->SetBinEdges(0,fBinsZv);
 610   fGenNDEventMatrix->SetBinEdges(1,fBinsMult);
 611   fGenNDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 612   fGenNDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 613   fGenNDEventMatrix->Sumw2();
 614
 615   fGenNSDEventMatrix = new THnSparseF("fGenNSDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 616   fGenNSDEventMatrix->SetBinEdges(0,fBinsZv);
 617   fGenNSDEventMatrix->SetBinEdges(1,fBinsMult);
 618   fGenNSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 619   fGenNSDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 620   fGenNSDEventMatrix->Sumw2();
 621
 622   //
 623   fTriggerEventMatrix = new THnSparseF("fTriggerEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 624   fTriggerEventMatrix->SetBinEdges(0,fBinsZv);
 625   fTriggerEventMatrix->SetBinEdges(1,fBinsMult);
 626   fTriggerEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 627   fTriggerEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 628   fTriggerEventMatrix->Sumw2();
 629
 630   fTriggerSDEventMatrix = new THnSparseF("fTriggerSDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 631   fTriggerSDEventMatrix->SetBinEdges(0,fBinsZv);
 632   fTriggerSDEventMatrix->SetBinEdges(1,fBinsMult);
 633   fTriggerSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 634   fTriggerSDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 635   fTriggerSDEventMatrix->Sumw2();
 636
 637   fTriggerDDEventMatrix = new THnSparseF("fTriggerDDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 638   fTriggerDDEventMatrix->SetBinEdges(0,fBinsZv);
 639   fTriggerDDEventMatrix->SetBinEdges(1,fBinsMult);
 640   fTriggerDDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 641   fTriggerDDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 642   fTriggerDDEventMatrix->Sumw2();
 643
 644   fTriggerNDEventMatrix = new THnSparseF("fTriggerNDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 645   fTriggerNDEventMatrix->SetBinEdges(0,fBinsZv);
 646   fTriggerNDEventMatrix->SetBinEdges(1,fBinsMult);
 647   fTriggerNDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 648   fTriggerNDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 649   fTriggerNDEventMatrix->Sumw2();
 650
 651   fTriggerNSDEventMatrix = new THnSparseF("fTriggerNSDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 652   fTriggerNSDEventMatrix->SetBinEdges(0,fBinsZv);
 653   fTriggerNSDEventMatrix->SetBinEdges(1,fBinsMult);
 654   fTriggerNSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 655   fTriggerNSDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 656   fTriggerNSDEventMatrix->Sumw2();
 657
 658   //
 659   fRecEventMatrix = new THnSparseF("fRecEventMatrix","mcZv:multTrueMC:Cent",3,binsEventMatrixCent); //,minEventMatrix,maxEventMatrix);
 660   fRecEventMatrix->SetBinEdges(0,fBinsZv);
 661   fRecEventMatrix->SetBinEdges(1,fBinsMult);
 662   fRecEventMatrix->SetBinEdges(2,fBinsCentrality);
 663   fRecEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 664   fRecEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 665   fRecEventMatrix->GetAxis(2)->SetTitle("Centrality");
 666   fRecEventMatrix->Sumw2();
 667
 668   fRecSDEventMatrix = new THnSparseF("fRecSDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 669   fRecSDEventMatrix->SetBinEdges(0,fBinsZv);
 670   fRecSDEventMatrix->SetBinEdges(1,fBinsMult);
 671   fRecSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 672   fRecSDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 673   fRecSDEventMatrix->Sumw2();
 674
 675   fRecDDEventMatrix = new THnSparseF("fRecDDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 676   fRecDDEventMatrix->SetBinEdges(0,fBinsZv);
 677   fRecDDEventMatrix->SetBinEdges(1,fBinsMult);
 678   fRecDDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 679   fRecDDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 680   fRecDDEventMatrix->Sumw2();
 681
 682   fRecNDEventMatrix = new THnSparseF("fRecNDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 683   fRecNDEventMatrix->SetBinEdges(0,fBinsZv);
 684   fRecNDEventMatrix->SetBinEdges(1,fBinsMult);
 685   fRecNDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 686   fRecNDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 687   fRecNDEventMatrix->Sumw2();
 688
 689   fRecNSDEventMatrix = new THnSparseF("fRecNSDEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 690   fRecNSDEventMatrix->SetBinEdges(0,fBinsZv);
 691   fRecNSDEventMatrix->SetBinEdges(1,fBinsMult);
 692   fRecNSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 693   fRecNSDEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 694   fRecNSDEventMatrix->Sumw2();
 695
 696   fRecCandleEventMatrix = new THnSparseF("fRecCandleEventMatrix","mcZv:multTrueMC",2,binsEventMatrix); //,minEventMatrix,maxEventMatrix);
 697   fRecCandleEventMatrix->SetBinEdges(0,fBinsZv);
 698   fRecCandleEventMatrix->SetBinEdges(1,fBinsMult);
 699   fRecCandleEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 700   fRecCandleEventMatrix->GetAxis(1)->SetTitle("true multiplicity (MC)");
 701   fRecCandleEventMatrix->Sumw2();
 702
 703   //
 704   // track to event corrections
 705   //
 706
 707   fGenTrackEventMatrix = new THnSparseF("fGenTrackEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 708   fGenTrackEventMatrix->SetBinEdges(0,fBinsZv);
 709   fGenTrackEventMatrix->SetBinEdges(1,fBinsPtCorr);
 710   fGenTrackEventMatrix->SetBinEdges(2,fBinsEta);
 711   fGenTrackEventMatrix->SetBinEdges(3,fBinsMultTE);
 712   fGenTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 713   fGenTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 714   fGenTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
 715   fGenTrackEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 716   fGenTrackEventMatrix->Sumw2();
 717
 718   fGenTrackSDEventMatrix = new THnSparseF("fGenTrackSDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 719   fGenTrackSDEventMatrix->SetBinEdges(0,fBinsZv);
 720   fGenTrackSDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 721   fGenTrackSDEventMatrix->SetBinEdges(2,fBinsEta);
 722   fGenTrackSDEventMatrix->SetBinEdges(3,fBinsMultTE);
 723   fGenTrackSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 724   fGenTrackSDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 725   fGenTrackSDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 726   fGenTrackSDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 727   fGenTrackSDEventMatrix->Sumw2();
 728
 729   fGenTrackDDEventMatrix = new THnSparseF("fGenTrackDDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 730   fGenTrackDDEventMatrix->SetBinEdges(0,fBinsZv);
 731   fGenTrackDDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 732   fGenTrackDDEventMatrix->SetBinEdges(2,fBinsEta);
 733   fGenTrackDDEventMatrix->SetBinEdges(3,fBinsMultTE);
 734   fGenTrackDDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 735   fGenTrackDDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 736   fGenTrackDDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 737   fGenTrackDDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 738   fGenTrackDDEventMatrix->Sumw2();
 739
 740   fGenTrackNDEventMatrix = new THnSparseF("fGenTrackNDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 741   fGenTrackNDEventMatrix->SetBinEdges(0,fBinsZv);
 742   fGenTrackNDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 743   fGenTrackNDEventMatrix->SetBinEdges(2,fBinsEta);
 744   fGenTrackNDEventMatrix->SetBinEdges(3,fBinsMultTE);
 745   fGenTrackNDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 746   fGenTrackNDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 747   fGenTrackNDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 748   fGenTrackNDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 749   fGenTrackNDEventMatrix->Sumw2();
 750
 751   fGenTrackNSDEventMatrix = new THnSparseF("fGenTrackNSDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 752   fGenTrackNSDEventMatrix->SetBinEdges(0,fBinsZv);
 753   fGenTrackNSDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 754   fGenTrackNSDEventMatrix->SetBinEdges(2,fBinsEta);
 755   fGenTrackNSDEventMatrix->SetBinEdges(3,fBinsMultTE);
 756   fGenTrackNSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 757   fGenTrackNSDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 758   fGenTrackNSDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 759   fGenTrackNSDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 760   fGenTrackNSDEventMatrix->Sumw2();
 761
 762
 763   //
 764   fTriggerTrackEventMatrix = new THnSparseF("fTriggerTrackEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 765   fTriggerTrackEventMatrix->SetBinEdges(0,fBinsZv);
 766   fTriggerTrackEventMatrix->SetBinEdges(1,fBinsPtCorr);
 767   fTriggerTrackEventMatrix->SetBinEdges(2,fBinsEta);
 768   fTriggerTrackEventMatrix->SetBinEdges(3,fBinsMultTE);
 769   fTriggerTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 770   fTriggerTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 771   fTriggerTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
 772   fTriggerTrackEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 773   fTriggerTrackEventMatrix->Sumw2();
 774
 775   fTriggerTrackSDEventMatrix = new THnSparseF("fTriggerTrackSDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 776   fTriggerTrackSDEventMatrix->SetBinEdges(0,fBinsZv);
 777   fTriggerTrackSDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 778   fTriggerTrackSDEventMatrix->SetBinEdges(2,fBinsEta);
 779   fTriggerTrackSDEventMatrix->SetBinEdges(3,fBinsMultTE);
 780   fTriggerTrackSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 781   fTriggerTrackSDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 782   fTriggerTrackSDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 783   fTriggerTrackSDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 784   fTriggerTrackSDEventMatrix->Sumw2();
 785
 786   fTriggerTrackDDEventMatrix = new THnSparseF("fTriggerTrackDDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 787   fTriggerTrackDDEventMatrix->SetBinEdges(0,fBinsZv);
 788   fTriggerTrackDDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 789   fTriggerTrackDDEventMatrix->SetBinEdges(2,fBinsEta);
 790   fTriggerTrackDDEventMatrix->SetBinEdges(3,fBinsMultTE);
 791   fTriggerTrackDDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 792   fTriggerTrackDDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 793   fTriggerTrackDDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 794   fTriggerTrackDDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 795   fTriggerTrackDDEventMatrix->Sumw2();
 796
 797   fTriggerTrackNDEventMatrix = new THnSparseF("fTriggerTrackNDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 798   fTriggerTrackNDEventMatrix->SetBinEdges(0,fBinsZv);
 799   fTriggerTrackNDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 800   fTriggerTrackNDEventMatrix->SetBinEdges(2,fBinsEta);
 801   fTriggerTrackNDEventMatrix->SetBinEdges(3,fBinsMultTE);
 802   fTriggerTrackNDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 803   fTriggerTrackNDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 804   fTriggerTrackNDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 805   fTriggerTrackNDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 806   fTriggerTrackNDEventMatrix->Sumw2();
 807
 808   fTriggerTrackNSDEventMatrix = new THnSparseF("fTriggerTrackNSDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 809   fTriggerTrackNSDEventMatrix->SetBinEdges(0,fBinsZv);
 810   fTriggerTrackNSDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 811   fTriggerTrackNSDEventMatrix->SetBinEdges(2,fBinsEta);
 812   fTriggerTrackNSDEventMatrix->SetBinEdges(3,fBinsMultTE);
 813   fTriggerTrackNSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 814   fTriggerTrackNSDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 815   fTriggerTrackNSDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 816   fTriggerTrackNSDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 817   fTriggerTrackNSDEventMatrix->Sumw2();
 818
 819
 820   //
 821   fRecTrackEventMatrix = new THnSparseF("fRecTrackEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 822   fRecTrackEventMatrix->SetBinEdges(0,fBinsZv);
 823   fRecTrackEventMatrix->SetBinEdges(1,fBinsPtCorr);
 824   fRecTrackEventMatrix->SetBinEdges(2,fBinsEta);
 825   fRecTrackEventMatrix->SetBinEdges(3,fBinsMultTE);
 826   fRecTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 827   fRecTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 828   fRecTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
 829   fRecTrackEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 830   fRecTrackEventMatrix->Sumw2();
 831
 832
 833   fRecTrackSDEventMatrix = new THnSparseF("fRecTrackSDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 834   fRecTrackSDEventMatrix->SetBinEdges(0,fBinsZv);
 835   fRecTrackSDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 836   fRecTrackSDEventMatrix->SetBinEdges(2,fBinsEta);
 837   fRecTrackSDEventMatrix->SetBinEdges(3,fBinsMultTE);
 838   fRecTrackSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 839   fRecTrackSDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 840   fRecTrackSDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 841   fRecTrackSDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 842   fRecTrackSDEventMatrix->Sumw2();
 843
 844   fRecTrackDDEventMatrix = new THnSparseF("fRecTrackDDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 845   fRecTrackDDEventMatrix->SetBinEdges(0,fBinsZv);
 846   fRecTrackDDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 847   fRecTrackDDEventMatrix->SetBinEdges(2,fBinsEta);
 848   fRecTrackDDEventMatrix->SetBinEdges(3,fBinsMultTE);
 849   fRecTrackDDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 850   fRecTrackDDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 851   fRecTrackDDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 852   fRecTrackDDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 853   fRecTrackDDEventMatrix->Sumw2();
 854
 855   fRecTrackNDEventMatrix = new THnSparseF("fRecTrackNDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 856   fRecTrackNDEventMatrix->SetBinEdges(0,fBinsZv);
 857   fRecTrackNDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 858   fRecTrackNDEventMatrix->SetBinEdges(2,fBinsEta);
 859   fRecTrackNDEventMatrix->SetBinEdges(3,fBinsMultTE);
 860   fRecTrackNDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 861   fRecTrackNDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 862   fRecTrackNDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 863   fRecTrackNDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 864   fRecTrackNDEventMatrix->Sumw2();
 865
 866   fRecTrackNSDEventMatrix = new THnSparseF("fRecTrackNSDEventMatrix","mcZv:mcPt:mcEta:multTrueMC",4,binsTrackEventCorrMatrix);
 867   fRecTrackNSDEventMatrix->SetBinEdges(0,fBinsZv);
 868   fRecTrackNSDEventMatrix->SetBinEdges(1,fBinsPtCorr);
 869   fRecTrackNSDEventMatrix->SetBinEdges(2,fBinsEta);
 870   fRecTrackNSDEventMatrix->SetBinEdges(3,fBinsMultTE);
 871   fRecTrackNSDEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 872   fRecTrackNSDEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 873   fRecTrackNSDEventMatrix->GetAxis(2)->SetTitle("mcEta");
 874   fRecTrackNSDEventMatrix->GetAxis(3)->SetTitle("true multiplicity (MC)");
 875   fRecTrackNSDEventMatrix->Sumw2();
 876   //
 877   // tracks correction matrices
 878   //
 879   //fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);
 880   fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);
 881   fGenTrackMatrix->SetBinEdges(0,fBinsZv);
 882   //fGenTrackMatrix->SetBinEdges(1,binsPt);
 883   fGenTrackMatrix->SetBinEdges(1,fBinsPtCorr);
 884   fGenTrackMatrix->SetBinEdges(2,fBinsEta);
 885   fGenTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 886   fGenTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 887   fGenTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 888   fGenTrackMatrix->Sumw2();
 889
 890   //fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);
 891   fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrixCent);
 892   fGenPrimTrackMatrix->SetBinEdges(0,fBinsZv);
 893   //fGenPrimTrackMatrix->SetBinEdges(1,binsPt);
 894   fGenPrimTrackMatrix->SetBinEdges(1,fBinsPtCorr);
 895   fGenPrimTrackMatrix->SetBinEdges(2,fBinsEta);
 896   fGenPrimTrackMatrix->SetBinEdges(3,fBinsCentrality);
 897   fGenPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 898   fGenPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 899   fGenPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 900   fGenPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 901   fGenPrimTrackMatrix->Sumw2();
 902
 903   //fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);
 904   fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrixCent);
 905   fRecPrimTrackMatrix->SetBinEdges(0,fBinsZv);
 906   //fRecPrimTrackMatrix->SetBinEdges(1,binsPt);
 907   fRecPrimTrackMatrix->SetBinEdges(1,fBinsPtCorr);
 908   fRecPrimTrackMatrix->SetBinEdges(2,fBinsEta);
 909   fRecPrimTrackMatrix->SetBinEdges(3,fBinsCentrality);
 910   fRecPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 911   fRecPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 912   fRecPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 913   fRecPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 914   fRecPrimTrackMatrix->Sumw2();
 915
 916   //
 917   //fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);
 918   fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrixCent);
 919   fRecTrackMatrix->SetBinEdges(0,fBinsZv);
 920   //fRecTrackMatrix->SetBinEdges(1,binsPt);
 921   fRecTrackMatrix->SetBinEdges(1,fBinsPtCorr);
 922   fRecTrackMatrix->SetBinEdges(2,fBinsEta);
 923   fRecTrackMatrix->SetBinEdges(3,fBinsCentrality);
 924   fRecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 925   fRecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)");
 926   fRecTrackMatrix->GetAxis(2)->SetTitle("Eta");
 927   fRecTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 928   fRecTrackMatrix->Sumw2();
 929
 930   //fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);
 931   fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrixCent);
 932   fRecSecTrackMatrix->SetBinEdges(0,fBinsZv);
 933   //fRecSecTrackMatrix->SetBinEdges(1,binsPt);
 934   fRecSecTrackMatrix->SetBinEdges(1,fBinsPtCorr);
 935   fRecSecTrackMatrix->SetBinEdges(2,fBinsEta);
 936   fRecSecTrackMatrix->SetBinEdges(3,fBinsCentrality);
 937   fRecSecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 938   fRecSecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)");
 939   fRecSecTrackMatrix->GetAxis(2)->SetTitle("Eta");
 940   fRecSecTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 941   fRecSecTrackMatrix->Sumw2();
 942
 943   //
 944   //fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);
 945   fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);
 946   fRecMultTrackMatrix->SetBinEdges(0,fBinsZv);
 947   //fRecMultTrackMatrix->SetBinEdges(1,binsPt);
 948   fRecMultTrackMatrix->SetBinEdges(1,fBinsPtCorr);
 949   fRecMultTrackMatrix->SetBinEdges(2,fBinsEta);
 950   fRecMultTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 951   fRecMultTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 952   fRecMultTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 953   fRecMultTrackMatrix->Sumw2();
 954
 955   //
 956   // Control analysis histograms
 957   //
 958
 959   Int_t binsMCEventHist1[3]={100,100,140};
 960   Double_t minMCEventHist1[3]={-0.1,-0.1,-35.};
 961   Double_t maxMCEventHist1[3]={0.1,0.1,35.};
 962   fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv",3,binsMCEventHist1,minMCEventHist1,maxMCEventHist1);
 963   fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)");
 964   fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)");
 965   fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)");
 966   fMCEventHist1->Sumw2();
 967
 968   //
 969   Int_t binsRecEventHist1[3]={100,100,140};
 970   Double_t minRecEventHist1[3]={-3.,-3.,-35.};
 971   Double_t maxRecEventHist1[3]={3.,3.,35.};
 972
 973   fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv",3,binsRecEventHist1,minRecEventHist1,maxRecEventHist1);
 974   fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)");
 975   fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)");
 976   fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)");
 977   fRecEventHist1->Sumw2();
 978
 979   //
 980   Int_t binsRecEventHist2[3]={fZvNbins,fMultNbins,fMultNbins};
 981 //   Double_t minRecEventHist2[3]={-25.,-0.5,-0.5};
 982 //   Double_t maxRecEventHist2[3]={25.,149.5,149.5};
 983
 984   fRecEventHist2 = new THnSparseF("fRecEventHist2","Zv:multMB:multRecMult",3,binsRecEventHist2); //,minRecEventHist2,maxRecEventHist2);
 985   fRecEventHist2->SetBinEdges(0,fBinsZv);
 986   fRecEventHist2->SetBinEdges(1,fBinsMult);
 987   fRecEventHist2->SetBinEdges(2,fBinsMult);
 988   fRecEventHist2->GetAxis(0)->SetTitle("Zv (cm)");
 989   fRecEventHist2->GetAxis(1)->SetTitle("multiplicity MB");
 990   fRecEventHist2->GetAxis(2)->SetTitle("multiplicity (multCuts)");
 991   fRecEventHist2->Sumw2();
 992
 993   //
 994   Double_t kFact = 0.1;
 995   Int_t binsRecMCEventHist1[3]={100,100,100};
 996   Double_t minRecMCEventHist1[3]={-10.0*kFact,-10.0*kFact,-10.0*kFact};
 997   Double_t maxRecMCEventHist1[3]={10.0*kFact,10.0*kFact,10.0*kFact};
 998
 999   fRecMCEventHist1 = new THnSparseF("fRecMCEventHist1","Xv-mcXv:Yv-mcYv:Zv-mcZv",3,binsRecMCEventHist1,minRecMCEventHist1,maxRecMCEventHist1);
1000   fRecMCEventHist1->GetAxis(0)->SetTitle("Xv-mcXv (cm)");
1001   fRecMCEventHist1->GetAxis(1)->SetTitle("Yv-mcYv (cm)");
1002   fRecMCEventHist1->GetAxis(2)->SetTitle("Zv-mcZv (cm)");
1003   fRecMCEventHist1->Sumw2();
1004
1005   //
1006   Int_t binsRecMCEventHist2[3]={100,100,fMultNbins};
1007   Double_t minRecMCEventHist2[3]={-10.0*kFact,-10.0*kFact,0.0};
1008   Double_t maxRecMCEventHist2[3]={10.0*kFact,10.0*kFact,149.50};
1009
1010   fRecMCEventHist2 = new THnSparseF("fRecMCEventHist2","Xv-mcXv:Zv-mcZv:multMB",3,binsRecMCEventHist2,minRecMCEventHist2,maxRecMCEventHist2);
1011   fRecMCEventHist2->SetBinEdges(2,fBinsMult);
1012   fRecMCEventHist2->GetAxis(0)->SetTitle("Xv-mcXv (cm)");
1013   fRecMCEventHist2->GetAxis(1)->SetTitle("Zv-mcZv (cm)");
1014   fRecMCEventHist2->GetAxis(2)->SetTitle("multiplicity MB");
1015   fRecMCEventHist2->Sumw2();
1016
1017   Int_t binsRecMCEventHist3[2]={fMultNbins,5};
1018   Double_t minRecMCEventHist3[2]={-0.5,0.0};
1019   Double_t maxRecMCEventHist3[2]={149.50,5.0};
1020   fRecMCEventHist3 = new THnSparseF("fRecMCEventHist3","multRecMult:EventType (ND, DD, SD)",2,binsRecMCEventHist3,minRecMCEventHist3,maxRecMCEventHist3);
1021   fRecMCEventHist3->SetBinEdges(0,fBinsMult);
1022   fRecMCEventHist3->GetAxis(0)->SetTitle("multiplicity (multCuts)");
1023   fRecMCEventHist3->GetAxis(1)->SetTitle("EventType");
1024   fRecMCEventHist3->Sumw2();
1025
1026   //
1027   char name[256];
1028   char title[256];
1029   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++)
1030   {
1031   // THnSparse track histograms
1032
1033   Int_t binsMCTrackHist1[3]=  {fPtCorrNbins, fEtaNbins, 90};
1034   Double_t minMCTrackHist1[3]={0.,-1.,0.};
1035   Double_t maxMCTrackHist1[3]={10.,1.,2.*TMath::Pi()};
1036   snprintf(name,256,"fMCTrackHist1_%d",i);
1037   snprintf(title,256,"mcPt:mcEta:mcPhi");
1038
1039   fMCTrackHist1[i] = new THnSparseF(name,title,3,binsMCTrackHist1,minMCTrackHist1,maxMCTrackHist1);
1040   fMCTrackHist1[i]->SetBinEdges(0,fBinsPtCorr);
1041   fMCTrackHist1[i]->SetBinEdges(1,fBinsEta);
1042   fMCTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
1043   fMCTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
1044   fMCTrackHist1[i]->GetAxis(2)->SetTitle("mcPhi (rad)");
1045   fMCTrackHist1[i]->Sumw2();
1046
1047   Int_t binsMCPrimTrackHist1[5]=  {fPtCorrNbins,fEtaNbins,6,20,4000};
1048   Double_t minMCPrimTrackHist1[5]={0.,-1.,0.,0.,0.};
1049   Double_t maxMCPrimTrackHist1[5]={10.,1.,6.,20.,4000.};
1050   snprintf(name,256,"fMCPrimTrackHist1_%d",i);
1051   snprintf(title,256,"mcPt:mcEta:pid:mech:mother");
1052
1053   fMCPrimTrackHist1[i] = new THnSparseF(name,title,5,binsMCPrimTrackHist1,minMCPrimTrackHist1,maxMCPrimTrackHist1);
1054   fMCPrimTrackHist1[i]->SetBinEdges(0,fBinsPtCorr);
1055   fMCPrimTrackHist1[i]->SetBinEdges(1,fBinsEta);
1056   fMCPrimTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
1057   fMCPrimTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
1058   fMCPrimTrackHist1[i]->GetAxis(2)->SetTitle("pid");
1059   fMCPrimTrackHist1[i]->GetAxis(3)->SetTitle("mech");
1060   fMCPrimTrackHist1[i]->GetAxis(4)->SetTitle("mother");
1061   fMCPrimTrackHist1[i]->Sumw2();
1062
1063   Int_t binsMCPrimTrackHist2[3]=  {4000,20,4000};
1064   Double_t minMCPrimTrackHist2[3]={0.,0.,0.};
1065   Double_t maxMCPrimTrackHist2[3]={4000.,20.,4000.};
1066   snprintf(name,256,"fMCPrimTrackHist2_%d",i);
1067   snprintf(title,256,"pdg:mech:mother");
1068
1069   fMCPrimTrackHist2[i] = new THnSparseF(name,title,3,binsMCPrimTrackHist2,minMCPrimTrackHist2,maxMCPrimTrackHist2);
1070   fMCPrimTrackHist2[i]->GetAxis(0)->SetTitle("pdg");
1071   fMCPrimTrackHist2[i]->GetAxis(1)->SetTitle("mech");
1072   fMCPrimTrackHist2[i]->GetAxis(2)->SetTitle("mother");
1073   fMCPrimTrackHist2[i]->Sumw2();
1074
1075   Int_t binsMCSecTrackHist1[5]=  {fPtCorrNbins,fEtaNbins,6,20,4000};
1076   Double_t minMCSecTrackHist1[5]={0.,-1.,0.,0.,0.};
1077   Double_t maxMCSecTrackHist1[5]={10.,1.,6.,20.,4000.};
1078   snprintf(name,256,"fMCSecTrackHist1_%d",i);
1079   snprintf(title,256,"mcPt:mcEta:pid:mech:mother");
1080
1081   fMCSecTrackHist1[i] = new THnSparseF(name,title,5,binsMCSecTrackHist1,minMCSecTrackHist1,maxMCSecTrackHist1);
1082   fMCSecTrackHist1[i]->SetBinEdges(0,fBinsPtCorr);
1083   fMCSecTrackHist1[i]->SetBinEdges(1,fBinsEta);
1084   fMCSecTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
1085   fMCSecTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
1086   fMCSecTrackHist1[i]->GetAxis(2)->SetTitle("pid");
1087   fMCSecTrackHist1[i]->GetAxis(3)->SetTitle("mech");
1088   fMCSecTrackHist1[i]->GetAxis(4)->SetTitle("mother");
1089   fMCSecTrackHist1[i]->Sumw2();
1090
1091   //
1092
1093   //
1094
1095   Int_t binsRecTrackHist1[3]={fPtNbins,fEtaNbins,90};
1096   Double_t minRecTrackHist1[3]={0.,-1.,0.};
1097   Double_t maxRecTrackHist1[3]={10.,1.,2.*TMath::Pi()};
1098   snprintf(name,256,"fRecTrackHist1_%d",i);
1099   snprintf(title,256,"Pt:Eta:Phi");
1100   fRecTrackHist1[i] = new THnSparseF(name,title,3,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1);
1101   fRecTrackHist1[i]->SetBinEdges(0,fBinsPt);
1102   fRecTrackHist1[i]->SetBinEdges(1,fBinsEta);
1103   fRecTrackHist1[i]->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
1104   fRecTrackHist1[i]->GetAxis(1)->SetTitle("#eta");
1105   fRecTrackHist1[i]->GetAxis(2)->SetTitle("#phi (rad)");
1106   fRecTrackHist1[i]->Sumw2();
1107
1108   //
1109   Int_t binsRecTrackMultHist1[2]={fPtNbins,fMultNbins};
1110   snprintf(name,256,"fRecTrackMultHist_%d",i);
1111   snprintf(title,256,"Pt:Mult");
1112   fRecTrackMultHist1[i] = new THnSparseF(name,title,2,binsRecTrackMultHist1); //,minRecTrackMultHist1,maxRecTrackMultHist1);
1113   fRecTrackMultHist1[i]->SetBinEdges(0,fBinsPt);
1114   fRecTrackMultHist1[i]->SetBinEdges(1,fBinsMult);
1115   fRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
1116   fRecTrackMultHist1[i]->GetAxis(1)->SetTitle("multiplicity");
1117   fRecTrackMultHist1[i]->Sumw2();
1118   }
1119
1120   Int_t binsRecMCTrackHist1[4] = {fPtCorrNbins,fEtaNbins,100,100};
1121   Double_t minRecMCTrackHist1[4]={0.,-1.,-0.5,-0.5};
1122   Double_t maxRecMCTrackHist1[4]={20.,1.,0.5,0.5};
1123   snprintf(name,256,"fRecMCTrackHist1");
1124   snprintf(title,256,"mcPt:mcEta:(Pt-mcPt)/mcPt:(Eta-mcEta)");
1125
1126   fRecMCTrackHist1 = new THnSparseF(name,title,4,binsRecMCTrackHist1,minRecMCTrackHist1,maxRecMCTrackHist1);
1127   fRecMCTrackHist1->SetBinEdges(0,fBinsPtCorr);
1128   fRecMCTrackHist1->SetBinEdges(1,fBinsEta);
1129   fRecMCTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
1130   fRecMCTrackHist1->GetAxis(1)->SetTitle("mcEta");
1131   fRecMCTrackHist1->GetAxis(2)->SetTitle("(Pt-mcPt)/mcPt");
1132   fRecMCTrackHist1->GetAxis(3)->SetTitle("Eta-mcEta");
1133
1134   Int_t binsMCMultRecTrackHist1[3] = {fPtCorrNbins,fEtaNbins,6};
1135   Double_t minMCMultRecTrackHist1[3]={0.,-1.,0.};
1136   Double_t maxMCMultRecTrackHist1[3]={20.,1.,6.};
1137   snprintf(name,256,"fMCMultRecTrackHist1");
1138   snprintf(title,256,"mcPt:mcEta:pid");
1139   fMCMultRecTrackHist1 = new THnSparseF(name,title,3,binsMCMultRecTrackHist1,minMCMultRecTrackHist1,maxMCMultRecTrackHist1);
1140   fMCMultRecTrackHist1->SetBinEdges(0,fBinsPtCorr);
1141   fMCMultRecTrackHist1->SetBinEdges(1,fBinsEta);
1142   fMCMultRecTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
1143   fMCMultRecTrackHist1->GetAxis(1)->SetTitle("mcEta");
1144   fMCMultRecTrackHist1->GetAxis(2)->SetTitle("pid");
1145
1146   //nClust:chi2PerClust:pt:eta:phi
1147   Int_t binsRecTrackHist2[5]={160,100,fPtNbins,fEtaNbins,90};
1148   Double_t minRecTrackHist2[5]={0., 0., 0, -1.5, 0.};
1149   Double_t maxRecRecTrackHist2[5]={160.,10., 16, 1.5, 2.*TMath::Pi()};
1150
1151   fRecTrackHist2 = new THnSparseF("fRecTrackHist2","nClust:chi2PerClust:pt:eta:phi",5,binsRecTrackHist2,minRecTrackHist2,maxRecRecTrackHist2);
1152   fRecTrackHist2->SetBinEdges(2,fBinsPt);
1153   fRecTrackHist2->SetBinEdges(3,fBinsEta);
1154   fRecTrackHist2->GetAxis(0)->SetTitle("nClust");
1155   fRecTrackHist2->GetAxis(1)->SetTitle("chi2PerClust");
1156   fRecTrackHist2->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
1157   fRecTrackHist2->GetAxis(3)->SetTitle("#eta");
1158   fRecTrackHist2->GetAxis(4)->SetTitle("#phi (rad)");
1159   fRecTrackHist2->Sumw2();
1160
1161   Int_t binsEventCount[3]={2,2,2};
1162   Double_t minEventCount[3]={0,0,0};
1163   Double_t maxEventCount[3]={2,2,2};
1164   fEventCount = new THnSparseF("fEventCount","trig vs trig+vertex",3,binsEventCount,minEventCount,maxEventCount);
1165   fEventCount->GetAxis(0)->SetTitle("trig");
1166   fEventCount->GetAxis(1)->SetTitle("trig+vert");
1167   fEventCount->GetAxis(2)->SetTitle("selected");
1168   fEventCount->Sumw2();
1169
1170   Int_t binsEventMultHist[3]={fMultNbins,fMultNbins,fMultNbins};
1171   fEventMultHist = new THnSparseF("fEventMultHist","multMB:multRecMult:multRec",3,binsEventMultHist);
1172   fEventMultHist->GetAxis(0)->SetTitle("multMB");
1173   fEventMultHist->GetAxis(1)->SetTitle("multRecMult");
1174   fEventMultHist->GetAxis(2)->SetTitle("multRec");
1175   fEventMultHist->SetBinEdges(0,fBinsMult);
1176   fEventMultHist->SetBinEdges(1,fBinsMult);
1177   fEventMultHist->SetBinEdges(2,fBinsMult);
1178   fEventMultHist->Sumw2();
1179
1180   //
1181   // create histograms for centrality
1182   // (only if fDimensionsCentralityEstimators > 0)
1183
1184 //   if (fDimensionsCentralityEstimators > 0) {
1185
1186     // fCentralityEventHist
1187     //zv:multRec:multMB:cent[n]
1188     Int_t dimsCentralityEventHist = fDimensionsCentralityEstimators + 3;
1189     Int_t* binsCentralityEventHist = new Int_t[dimsCentralityEventHist];
1190     TString titleCentralityEventHist("zV:multRecMult:multMB");
1191     binsCentralityEventHist[0] = fZvNbins;
1192     binsCentralityEventHist[1] = fMultNbins;
1193     binsCentralityEventHist[2] = fMultNbins;
1194     for (Int_t i=1; i<=fDimensionsCentralityEstimators ; i++) {
1195       binsCentralityEventHist[i+2] = fCentralityNbins;
1196       titleCentralityEventHist += ":cent";
1197       titleCentralityEventHist += TString(GetCentralityEstimator(i));
1198     }
1199     fCentralityEventHist = new THnSparseF("fCentralityEventHist",titleCentralityEventHist.Data(),dimsCentralityEventHist,binsCentralityEventHist);
1200     fCentralityEventHist->SetBinEdges(0,fBinsZv);
1201     fCentralityEventHist->SetBinEdges(1,fBinsMult);
1202     fCentralityEventHist->SetBinEdges(2,fBinsMult);
1203     fCentralityEventHist->GetAxis(0)->SetTitle("Zv (cm)");
1204     fCentralityEventHist->GetAxis(1)->SetTitle("multRecMult");
1205     fCentralityEventHist->GetAxis(2)->SetTitle("multMB");
1206     for (Int_t i=1; i<=fDimensionsCentralityEstimators ; i++) {
1207       fCentralityEventHist->SetBinEdges(i+2,fBinsCentrality);
1208       TString axistitle(GetCentralityEstimator(i));
1209       axistitle += " Centrality";
1210       fCentralityEventHist->GetAxis(i+2)->SetTitle(axistitle);
1211     }
1212     fCentralityEventHist->Sumw2();
1213
1214     // fCentralityTrackHist
1215     //zv:pt:eta:multRec:multMB:cent[n]
1216     Int_t dimsCentralityTrackHist = fDimensionsCentralityEstimators + 5;
1217     Int_t* binsCentralityTrackHist = new Int_t[dimsCentralityTrackHist];
1218     TString titleCentralityTrackHist("zV:pT:eta:multRecMult:multMB");
1219     binsCentralityTrackHist[0] = fZvNbins;
1220     binsCentralityTrackHist[1] = fPtNbins;
1221     binsCentralityTrackHist[2] = fEtaNbins;
1222     binsCentralityTrackHist[3] = fMultNbins;
1223     binsCentralityTrackHist[4] = fMultNbins;
1224     for (Int_t i=1; i<=fDimensionsCentralityEstimators ; i++) {
1225       binsCentralityTrackHist[i+4] = fCentralityNbins;
1226       titleCentralityTrackHist += ":cent";
1227       titleCentralityTrackHist += TString(GetCentralityEstimator(i));
1228     }
1229     fCentralityTrackHist = new THnSparseF("fCentralityTrackHist",titleCentralityTrackHist.Data(),dimsCentralityTrackHist,binsCentralityTrackHist);
1230     fCentralityTrackHist->SetBinEdges(0,fBinsZv);
1231     fCentralityTrackHist->SetBinEdges(1,fBinsPt);
1232     fCentralityTrackHist->SetBinEdges(2,fBinsEta);
1233     fCentralityTrackHist->SetBinEdges(3,fBinsMult);
1234     fCentralityTrackHist->SetBinEdges(4,fBinsMult);
1235     fCentralityTrackHist->GetAxis(0)->SetTitle("Zv (cm)");
1236     fCentralityTrackHist->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
1237     fCentralityTrackHist->GetAxis(2)->SetTitle("#eta");
1238     fCentralityTrackHist->GetAxis(3)->SetTitle("multRecMult");
1239     fCentralityTrackHist->GetAxis(4)->SetTitle("multMB");
1240     for (Int_t i=1; i<=fDimensionsCentralityEstimators ; i++) {
1241       fCentralityTrackHist->SetBinEdges(i+4,fBinsCentrality);
1242       TString axistitle(GetCentralityEstimator(i));
1243       axistitle += " Centrality";
1244       fCentralityTrackHist->GetAxis(i+4)->SetTitle(axistitle);
1245     }
1246     fCentralityTrackHist->Sumw2();
1247
1248     fNVCentralityEvent = fDimensionsCentralityEstimators + 3;
1249     fNVCentralityTrack = fDimensionsCentralityEstimators + 5;
1250     fVCentralityEvent = new Double_t[fNVCentralityEvent];
1251     fVCentralityTrack = new Double_t[fNVCentralityTrack];
1252  // }
1253
1254   // init folder
1255   fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");
1256
1257   fUtils = new AliAnalysisUtils();
1258
1259   // set init flag
1260   fIsInit = kTRUE;
1261 }
1262
1263 //_____________________________________________________________________________
1264 void AlidNdPtAnalysispPb::Process(AliESDEvent *const esdEvent, AliMCEvent *const mcEvent)
1265 {
1266   //  init if not done already
1267   if (!fIsInit) { Init(); }
1268
1269
1270   //
1271   // Process real and/or simulated events
1272   //
1273   if(!esdEvent) {
1274     AliDebug(AliLog::kError, "esdEvent not available");
1275     return;
1276   }
1277
1278   // get selection cuts
1279   AlidNdPtEventCuts *evtCuts = GetEventCuts();
1280   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts();
1281   AliESDtrackCuts *esdTrackCuts = GetTrackCuts();
1282   AlidNdPtAcceptanceCuts *multAccCuts = GetMultAcceptanceCuts();
1283
1284   if(!evtCuts || !accCuts  || !esdTrackCuts) {
1285     AliDebug(AliLog::kError, "cuts not available");
1286     return;
1287   }
1288
1289   // trigger selection
1290   Bool_t isEventTriggered = kTRUE;
1291   AliPhysicsSelection *physicsSelection = NULL;
1292   AliTriggerAnalysis* triggerAnalysis = NULL;
1293
1294   //
1295   AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
1296   if (!inputHandler)
1297   {
1298     Printf("ERROR: Could not receive input handler");
1299     return;
1300   }
1301
1302   if(evtCuts->IsTriggerRequired())
1303   {
1304     // always MB
1305     //isEventTriggered = inputHandler->IsEventSelected() & AliVEvent::kMB;
1306     isEventTriggered = inputHandler->IsEventSelected() & GetTriggerMask();
1307
1308     physicsSelection = static_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());
1309     if(!physicsSelection) return;
1310     //SetPhysicsTriggerSelection(physicsSelection);
1311
1312     if (isEventTriggered && (GetTrigger() == AliTriggerAnalysis::kV0AND)) {
1313       // set trigger (V0AND)
1314       triggerAnalysis = physicsSelection->GetTriggerAnalysis();
1315       if(!triggerAnalysis) return;
1316       isEventTriggered = triggerAnalysis->IsOfflineTriggerFired(esdEvent, GetTrigger());
1317     }
1318
1319     // calculate LHC background
1320     if(!IsUseMCInfo())
1321     {
1322       //
1323       // 0-multiplicity bin for LHC background correction
1324       //
1325       /* bin0 done in the train
1326       if( GetAnalysisMode() == AlidNdPtHelper::kTPCITS ||
1327           GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtx ||
1328           GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtxUpdate ||
1329           GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtx ||
1330           GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt )
1331       {
1332          physicsSelection->SetBin0CallbackViaPointer(&AlidNdPtAnalysispPb::IsBinZeroTrackSPDvtx);
1333       } else {
1334          physicsSelection->SetBin0CallbackViaPointer(&AlidNdPtAnalysispPb::IsBinZeroSPDvtx);
1335       }
1336       */
1337     }
1338   }
1339
1340
1341   // use MC information
1342   AliHeader* header = 0;
1343   AliGenEventHeader* genHeader = 0;
1344   AliStack* stack = 0;
1345   TArrayF vtxMC(3);
1346   AliPWG0Helper::MCProcessType evtType = AliPWG0Helper::kInvalidProcess;
1347
1348   Int_t multMCTrueTracks = 0;
1349   if(IsUseMCInfo())
1350   {
1351     //
1352     if(!mcEvent) {
1353       AliDebug(AliLog::kError, "mcEvent not available");
1354       return;
1355     }
1356     // get MC event header
1357     header = mcEvent->Header();
1358     if (!header) {
1359       AliDebug(AliLog::kError, "Header not available");
1360       return;
1361     }
1362     // MC particle stack
1363     stack = mcEvent->Stack();
1364     if (!stack) {
1365       AliDebug(AliLog::kError, "Stack not available");
1366       return;
1367     }
1368     // get event type (ND=0x1, DD=0x2, SD=0x4)
1369     evtType = AlidNdPtHelper::GetEventProcessTypePA(header);
1370     AliDebug(AliLog::kDebug+1, Form("Found process type %d", evtType));
1371
1372     // get MC vertex
1373     genHeader = header->GenEventHeader();
1374     if (!genHeader) {
1375       AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
1376       return;
1377     }
1378     genHeader->PrimaryVertex(vtxMC);
1379
1380     Double_t vMCEventHist1[3]={vtxMC[0],vtxMC[1],vtxMC[2]};
1381     fMCEventHist1->Fill(vMCEventHist1);
1382
1383     // multipliticy of all MC primary tracks
1384     // in Zv, pt and eta ranges)
1385     if (fRapidityShift == 0) {
1386       multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,multAccCuts);
1387     } else {
1388       multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,multAccCuts,fRapidityShift);
1389     }
1390
1391   } // end bUseMC
1392
1393
1394   // get reconstructed vertex
1395   const AliESDVertex* vtxESD = 0;
1396   Bool_t isRecVertex = kFALSE;
1397   /* the following lines have been commented out for same vertex as dNdEta analysis
1398   if(evtCuts->IsRecVertexRequired())
1399   {
1400     Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();
1401     Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();
1402     vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints);
1403     if(!vtxESD) return;
1404     isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, esdEvent->GetPrimaryVertexSPD(), GetAnalysisMode(), kFALSE);
1405   }
1406
1407   if( IsUseMCInfo() && !evtCuts->IsRecVertexRequired() ) {
1408     vtxESD = new AliESDVertex(vtxMC[2],10.,genHeader->NProduced(),"smearMC");
1409     if(!vtxESD) return;
1410     isRecVertex = kTRUE;
1411   }
1412   if(!vtxESD) return;
1413   */
1414   // the following lines have  been added for same vertex as dNdEta analysis
1415     vtxESD = esdEvent->GetPrimaryVertexTracks();
1416     if(!vtxESD || (vtxESD->GetNContributors()<1)) {
1417       // SPD vertex
1418       vtxESD = esdEvent->GetPrimaryVertexSPD();
1419       if (vtxESD->GetNContributors()>0) {
1420       TString vtxTyp = vtxESD->GetTitle();
1421       if ( !vtxTyp.Contains("vertexer: Z") || (vtxESD->GetDispersion()<0.04 && vtxESD->GetZRes()<0.25)) isRecVertex = kTRUE;
1422       }
1423
1424     } else {
1425        isRecVertex = kTRUE;
1426     }
1427
1428   // the previous lines have been added for same vertex as dNdEta analysis
1429
1430   Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex;
1431
1432   Bool_t isEventSelected2013 = kTRUE;
1433
1434   // selection and pileup rejection for 2013 p-A
1435   if (fIs2013pA) {
1436        if (fUtils->IsFirstEventInChunk(esdEvent)) { isEventSelected2013 = kFALSE;  }
1437        if (!fUtils->IsVertexSelected2013pA(esdEvent)) { isEventSelected2013 = kFALSE;  }
1438        if (fUtils->IsPileUpEvent(esdEvent)) { isEventSelected2013 = kFALSE;  }
1439   }
1440   isEventOK = isEventOK && isEventSelected2013;
1441   //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);
1442   //printf("GetAnalysisMode() %d \n",GetAnalysisMode());
1443
1444   Bool_t isTrigAndVertex = isEventTriggered && isEventOK;
1445
1446   Double_t vEventCount[3] = { (isEventTriggered && kTRUE) , isTrigAndVertex,  isTrigAndVertex && (TMath::Abs(vtxESD->GetZv()) < 10.)};
1447   fEventCount->Fill(vEventCount);
1448
1449   // vertex contributors
1450   Int_t multMBTracks = 0;
1451   if(GetAnalysisMode() == AlidNdPtHelper::kTPC)
1452   {
1453      if(vtxESD && vtxESD->GetStatus() && isRecVertex)
1454        multMBTracks = AlidNdPtHelper::GetTPCMBTrackMult(esdEvent,evtCuts,accCuts,esdTrackCuts);
1455   }
1456   else if( GetAnalysisMode() == AlidNdPtHelper::kTPCSPDvtx ||
1457            GetAnalysisMode()==AlidNdPtHelper::kTPCSPDvtxUpdate )
1458   {
1459      const AliMultiplicity* mult = esdEvent->GetMultiplicity();
1460      if(mult && vtxESD && vtxESD->GetStatus() && isRecVertex)
1461        multMBTracks = mult->GetNumberOfTracklets();
1462   }
1463   else if( GetAnalysisMode() == AlidNdPtHelper::kTPCITS ||
1464            GetAnalysisMode()==AlidNdPtHelper::kTPCITSHybrid ||
1465            GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtx ||
1466            GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtxUpdate ||
1467            GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtx ||
1468            GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt ||
1469            GetAnalysisMode() == AlidNdPtHelper::kITSStandAloneTrackSPDvtx ||
1470            GetAnalysisMode() == AlidNdPtHelper::kITSStandAloneTPCTrackSPDvtx
1471            )
1472   {
1473      if(vtxESD && vtxESD->GetStatus() && isRecVertex)
1474        multMBTracks = vtxESD->GetNContributors();
1475   }
1476   else {
1477     AliDebug(AliLog::kError, Form("Found analysis type %d", GetAnalysisMode()));
1478     return;
1479   }
1480
1481   Bool_t isEventSelected = kTRUE;
1482   if(evtCuts->IsEventSelectedRequired())
1483   {
1484      // select events with at least
1485      // one prompt track in acceptance
1486      // pT>0.5 GeV/c, |eta|<0.8 for the Cross Section studies
1487
1488      isEventSelected = AlidNdPtHelper::SelectEvent(esdEvent,esdTrackCuts);
1489      //printf("isEventSelected %d \n", isEventSelected);
1490   }
1491
1492   TObjArray *allChargedTracks=0;
1493   //Int_t multAll=0, multAcc=0, multRec=0;
1494   Int_t multAll=0, multRec=0;
1495   Int_t *labelsAll=0, *labelsAcc=0, *labelsRec=0;
1496
1497
1498   // check event cuts
1499   Int_t multRecMult=0;
1500   Double_t centralityD = -1;
1501   if(isEventOK && isEventTriggered && isEventSelected)
1502   {
1503     // get all charged tracks
1504     allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());
1505     if(!allChargedTracks) return;
1506
1507     Int_t entries = allChargedTracks->GetEntries();
1508     //printf("entries %d \n",entries);
1509
1510
1511     // calculate mult of reconstructed tracks
1512     for(Int_t i=0; i<entries;++i)
1513     {
1514       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
1515       if(!track) continue;
1516       if(track->Charge()==0) continue;
1517
1518
1519       // only postive charged
1520       if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0)
1521         continue;
1522
1523       // only negative charged
1524       if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0)
1525         continue;
1526
1527       if(GetMultTrackCuts()->AcceptTrack(track))
1528       {
1529           if(GetMultAcceptanceCuts()->AcceptTrack(track)) { multRecMult++; }
1530       }
1531     }
1532
1533
1534      // fill fCentralityEventHist only if necessary
1535      // also filles event-related value for fCentralityTrackHist
1536      if (fDimensionsCentralityEstimators > 0) {
1537        fVCentralityEvent[0] = vtxESD->GetZv();
1538        fVCentralityEvent[1] = multRecMult;
1539        fVCentralityEvent[2] = multMBTracks;
1540        fVCentralityTrack[0] = vtxESD->GetZv();
1541        fVCentralityTrack[3] = multRecMult;
1542        fVCentralityTrack[4] = multMBTracks;
1543        for (Int_t i=1; i<=fDimensionsCentralityEstimators ; i++) {
1544          // centrality determination
1545          Float_t centralityF = -1.;
1546          AliCentrality *esdCentrality = esdEvent->GetCentrality();
1547          centralityF = esdCentrality->GetCentralityPercentile(GetCentralityEstimator(i));
1548          fVCentralityEvent[i+2] = centralityF;
1549          fVCentralityTrack[i+4] = centralityF;
1550          if (1==i) { centralityD = centralityF; }
1551        }
1552        fCentralityEventHist->Fill(fVCentralityEvent);
1553      }
1554
1555
1556     labelsAll = new Int_t[entries];
1557     labelsAcc = new Int_t[entries];
1558     labelsRec = new Int_t[entries];
1559     for(Int_t i=0; i<entries;++i)
1560     {
1561       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
1562       if(!track) continue;
1563       if(track->Charge()==0) continue;
1564
1565       // only postive charged
1566       if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0)
1567         continue;
1568
1569       // only negative charged
1570       if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0)
1571         continue;
1572
1573       FillHistograms(track,stack,vtxESD->GetZv(),AlidNdPtHelper::kAllTracks, multRecMult);
1574       labelsAll[multAll] = TMath::Abs(track->GetLabel());
1575       multAll++;
1576
1577       //if(accCuts->AcceptTrack(track)) {
1578       //FillHistograms(track,stack,AlidNdPtHelper::kAccTracks);
1579       //labelsAcc[multAcc] = TMath::Abs(track->GetLabel());
1580       //multAcc++;
1581       //}
1582
1583       // esd track selection
1584       // ITS stand alone
1585       if(GetAnalysisMode() == AlidNdPtHelper::kITSStandAloneTrackSPDvtx)
1586       {
1587         if(!(track->GetStatus() & AliESDtrack::kITSpureSA)) continue;
1588         if(!(track->GetStatus() & AliESDtrack::kITSrefit)) continue;
1589         if(track->GetNcls(0)<4) continue;
1590         if(!track->HasPointOnITSLayer(0) && !track->HasPointOnITSLayer(1)) continue;
1591       }
1592       else if(GetAnalysisMode() == AlidNdPtHelper::kITSStandAloneTPCTrackSPDvtx)
1593       {
1594         //
1595         // ITS and TPC stand alone tracks
1596         //
1597         if(!(track->GetStatus() & AliESDtrack::kITSpureSA)) continue;
1598         if(!(track->GetStatus() & AliESDtrack::kITSrefit)) continue;
1599         if(track->GetNcls(0)<4) continue;
1600         if(!track->HasPointOnITSLayer(0) && !track->HasPointOnITSLayer(1)) continue;
1601
1602         // Check matching with TPC only track
1603         Bool_t hasMatch = kFALSE;
1604         for(Int_t j=0; j<entries;++j)
1605         {
1606           if(i==j) continue;
1607           AliESDtrack *track2 = (AliESDtrack*)allChargedTracks->At(j);
1608           if(!track2) continue;
1609           if(track2->Charge()==0) continue;
1610           if (!track2->GetTPCInnerParam()) continue;
1611
1612           // check loose cuts for TPC tracks
1613           if(!esdTrackCuts->AcceptTrack(track2)) continue;
1614
1615           AliExternalTrackParam *innerTPC = new AliExternalTrackParam(*(track2->GetTPCInnerParam()));
1616           if(!innerTPC) continue;
1617           Double_t x2[3]; track2->GetXYZ(x2);
1618           Double_t b2[3]; AliTracker::GetBxByBz(x2,b2);
1619           Double_t dz[2],cov[3];
1620           Bool_t isPropOK = innerTPC->PropagateToDCABxByBz(vtxESD,b2,kVeryBig,dz,cov);
1621           if(!isPropOK && innerTPC) {delete innerTPC; continue;}
1622
1623           // check matching
1624           if (TMath::Abs(track->GetY() - innerTPC->GetY()) > 3) { delete innerTPC; continue; }
1625           if (TMath::Abs(track->GetSnp() - innerTPC->GetSnp()) > 0.2) { delete innerTPC; continue; }
1626           if (TMath::Abs(track->GetTgl() - innerTPC->GetTgl()) > 0.2) { delete innerTPC; continue; }
1627
1628           hasMatch = kTRUE;
1629           if(innerTPC) delete innerTPC;
1630         }
1631
1632         if(!hasMatch) continue;
1633       }
1634       else {
1635         // check track cuts
1636
1637         if(!esdTrackCuts->AcceptTrack(track))
1638           continue;
1639       }
1640
1641       //
1642       Bool_t isOK = kFALSE;
1643       Double_t x[3]; track->GetXYZ(x);
1644       Double_t b[3]; AliTracker::GetBxByBz(x,b);
1645
1646       //
1647       // if TPC-ITS hybrid tracking (kTPCITSHybrid)
1648       // replace track parameters with TPC-ony track parameters
1649       //
1650       if( GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybrid ||
1651           GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtx ||
1652           GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt)
1653       {
1654         // Relate TPC-only tracks to Track or SPD vertex
1655         isOK = track->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig);
1656         if(!isOK) continue;
1657
1658         // replace esd track parameters with TPCinner
1659         AliExternalTrackParam  *tpcTrack  = new AliExternalTrackParam(*(track->GetTPCInnerParam()));
1660         if (tpcTrack) {
1661           track->Set(tpcTrack->GetX(),tpcTrack->GetAlpha(),tpcTrack->GetParameter(),tpcTrack->GetCovariance());
1662         }
1663         if(tpcTrack) delete tpcTrack;
1664       }
1665
1666
1667          // update track parameters using vertex point
1668          if( GetAnalysisMode() == AlidNdPtHelper::kTPCSPDvtxUpdate ||
1669              GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtxUpdate )
1670          {
1671            // update track parameters
1672              AliExternalTrackParam cParam;
1673              isOK = track->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig, &cParam);
1674              if(!isOK) continue;
1675              track->Set(cParam.GetX(),cParam.GetAlpha(),cParam.GetParameter(),cParam.GetCovariance());
1676
1677              if(accCuts->AcceptTrack(track)) {
1678                FillHistograms(track,stack,vtxESD->GetZv(),AlidNdPtHelper::kRecTracks,multRecMult);
1679                labelsRec[multRec] = TMath::Abs(track->GetLabel());
1680                multRec++;
1681              }
1682           }
1683           else {
1684              if(accCuts->AcceptTrack(track))
1685              {
1686                FillHistograms(track,stack,vtxESD->GetZv(),AlidNdPtHelper::kRecTracks,multRecMult);
1687                labelsRec[multRec] = TMath::Abs(track->GetLabel());
1688                multRec++;
1689              }
1690           }
1691      }
1692
1693
1694      // fill track multiplicity histograms
1695      // terribly slow
1696      // FillHistograms(allChargedTracks,labelsAll,multAll,labelsAcc,multAcc,labelsRec,multRec);
1697
1698      Double_t vRecEventHist1[3] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv()};
1699      fRecEventHist1->Fill(vRecEventHist1);
1700
1701      Double_t vRecEventHist2[3] = {vtxESD->GetZv(),multMBTracks,multRecMult};
1702      fRecEventHist2->Fill(vRecEventHist2);
1703
1704      //
1705      //Double_t vRecEventHist[2] = {vtxESD->GetZv(),multMBTracks};
1706      Double_t vRecEventHist[3] = {vtxESD->GetZv(),multMBTracks,centralityD};
1707      fRecEventHist->Fill(vRecEventHist);
1708
1709      // fill fEventMultHist for cross checks
1710      Double_t vEventMultHist[3] = {multMBTracks,multRecMult,multRec};
1711      fEventMultHist->Fill(vEventMultHist);
1712    }
1713
1714    if(IsUseMCInfo())
1715    {
1716      if(mcEvent) {
1717
1718      if(evtCuts->IsEventSelectedRequired())
1719      {
1720        // select events with at least
1721        // one MC primary track in acceptance
1722        // pT>0.5 GeV/c, |eta|<0.8 for the Cross Section studies
1723
1724        Bool_t isMCEventSelected = AlidNdPtHelper::SelectMCEvent(mcEvent);
1725        //printf("isMCEventSelected %d \n", isMCEventSelected);
1726        if(!isMCEventSelected) {
1727
1728         if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;
1729         if(labelsAll) delete [] labelsAll; labelsAll = 0;
1730         if(labelsAcc) delete [] labelsAcc; labelsAcc = 0;
1731         if(labelsRec) delete [] labelsRec; labelsRec = 0;
1732
1733         return;
1734        }
1735      }
1736
1737      Double_t vMultTrueEventMatrix[3] = { multRecMult, multMCTrueTracks, multMBTracks};
1738      if(isEventOK && isEventTriggered) {
1739        if(TMath::Abs(vtxMC[2]) < 10.0) // both Rec. and corresponding MC events must be accepted
1740          fEventMultCorrelationMatrix->Fill(vMultTrueEventMatrix);
1741      }
1742
1743      //
1744      // event level corrections (zv,N_MB)
1745      //
1746
1747      // all inelastic
1748      //Double_t vEventMatrix[2] = {vtxMC[2],multMBTracks};
1749      Double_t vEventMatrix[3] = {vtxMC[2],multMCTrueTracks,centralityD};
1750      fGenEventMatrix->Fill(vEventMatrix);
1751      if(isEventTriggered) fTriggerEventMatrix->Fill(vEventMatrix);
1752      if(isEventOK && isEventTriggered) fRecEventMatrix->Fill(vEventMatrix);
1753
1754      // single diffractive
1755      if(evtType == AliPWG0Helper::kSD) {
1756        fGenSDEventMatrix->Fill(vEventMatrix);
1757        if(isEventTriggered) fTriggerSDEventMatrix->Fill(vEventMatrix);
1758        if(isEventOK && isEventTriggered) fRecSDEventMatrix->Fill(vEventMatrix);
1759      }
1760
1761      // double diffractive
1762      if(evtType == AliPWG0Helper::kDD) {
1763        fGenDDEventMatrix->Fill(vEventMatrix);
1764        if(isEventTriggered) fTriggerDDEventMatrix->Fill(vEventMatrix);
1765        if(isEventOK && isEventTriggered)  fRecDDEventMatrix->Fill(vEventMatrix);
1766      }
1767
1768      // non diffractive
1769      if(evtType == AliPWG0Helper::kND) {
1770        fGenNDEventMatrix->Fill(vEventMatrix);
1771        if(isEventTriggered) fTriggerNDEventMatrix->Fill(vEventMatrix);
1772        if(isEventOK && isEventTriggered) fRecNDEventMatrix->Fill(vEventMatrix);
1773      }
1774
1775      // non single diffractive
1776      if(evtType != AliPWG0Helper::kSD) {
1777        fGenNSDEventMatrix->Fill(vEventMatrix);
1778        if(isEventTriggered) fTriggerNSDEventMatrix->Fill(vEventMatrix);
1779        if(isEventOK && isEventTriggered) fRecNSDEventMatrix->Fill(vEventMatrix);
1780      }
1781
1782      //
1783      // track-event level corrections (zv,pt,eta)
1784      //
1785      for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
1786      {
1787        TParticle* particle = stack->Particle(iMc);
1788        if (!particle)
1789        continue;
1790
1791        // only charged particles
1792        if(!particle->GetPDG()) continue;
1793        Double_t charge = particle->GetPDG()->Charge()/3.;
1794        if ( TMath::Abs(charge) < 0.001 )
1795         continue;
1796
1797        // only postive charged
1798        if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.)
1799         continue;
1800
1801        // only negative charged
1802        if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.)
1803        continue;
1804
1805        // physical primary
1806        Bool_t prim = stack->IsPhysicalPrimary(iMc);
1807        if(!prim) continue;
1808
1809        // checked accepted
1810        if(accCuts->AcceptTrack(particle))
1811        {
1812          //Double_t vTrackEventMatrix[3] = {vtxMC[2], particle->Pt(), particle->Eta()};
1813          Double_t vTrackEventMatrix[4] = {vtxMC[2],particle->Pt(),particle->Eta(),multMCTrueTracks};
1814          fGenTrackEventMatrix->Fill(vTrackEventMatrix);
1815
1816          if(evtType == AliPWG0Helper::kSD) {
1817            fGenTrackSDEventMatrix->Fill(vTrackEventMatrix);
1818          }
1819          if(evtType == AliPWG0Helper::kDD) {
1820            fGenTrackDDEventMatrix->Fill(vTrackEventMatrix);
1821          }
1822          if(evtType == AliPWG0Helper::kND) {
1823            fGenTrackNDEventMatrix->Fill(vTrackEventMatrix);
1824          }
1825          if(evtType != AliPWG0Helper::kSD) {
1826            fGenTrackNSDEventMatrix->Fill(vTrackEventMatrix);
1827          }
1828
1829          //
1830          if(!isEventTriggered) continue;
1831
1832          fTriggerTrackEventMatrix->Fill(vTrackEventMatrix);
1833          if(evtType == AliPWG0Helper::kSD) {
1834            fTriggerTrackSDEventMatrix->Fill(vTrackEventMatrix);
1835          }
1836          if(evtType == AliPWG0Helper::kDD) {
1837            fTriggerTrackDDEventMatrix->Fill(vTrackEventMatrix);
1838          }
1839          if(evtType == AliPWG0Helper::kND) {
1840            fTriggerTrackNDEventMatrix->Fill(vTrackEventMatrix);
1841          }
1842          if(evtType != AliPWG0Helper::kSD) {
1843            fTriggerTrackNSDEventMatrix->Fill(vTrackEventMatrix);
1844          }
1845
1846          //
1847          if(!isEventOK) continue;
1848
1849          fRecTrackEventMatrix->Fill(vTrackEventMatrix);
1850          if(evtType == AliPWG0Helper::kSD) {
1851            fRecTrackSDEventMatrix->Fill(vTrackEventMatrix);
1852          }
1853          if(evtType == AliPWG0Helper::kDD) {
1854            fRecTrackDDEventMatrix->Fill(vTrackEventMatrix);
1855          }
1856          if(evtType == AliPWG0Helper::kND) {
1857            fRecTrackNDEventMatrix->Fill(vTrackEventMatrix);
1858          }
1859          if(evtType != AliPWG0Helper::kSD) {
1860            fRecTrackNSDEventMatrix->Fill(vTrackEventMatrix);
1861          }
1862        }
1863      }
1864
1865      //
1866      // track-level corrections (zv,pt,eta)
1867      //
1868      if(isEventOK && isEventTriggered)
1869      {
1870
1871        // fill MC and rec event control histograms
1872        if(fHistogramsOn) {
1873          Double_t vRecMCEventHist1[3] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2]};
1874          fRecMCEventHist1->Fill(vRecMCEventHist1);
1875
1876          Double_t vRecMCEventHist2[3] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetZv()-vtxMC[2],multMBTracks};
1877          fRecMCEventHist2->Fill(vRecMCEventHist2);
1878
1879          Double_t vRecMCEventHist3[2] = {multRec,evtType};
1880          fRecMCEventHist3->Fill(vRecMCEventHist3);
1881        }
1882
1883        //
1884        // MC histograms for track efficiency studies
1885        //
1886        Int_t countRecCandle = 0;
1887        for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
1888        {
1889          TParticle* particle = stack->Particle(iMc);
1890          if (!particle)
1891          continue;
1892
1893          //Double_t vTrackMatrix[3] = {vtxMC[2],particle->Pt(),particle->Eta()}; .
1894          Double_t vTrackMatrix[4] = {vtxMC[2],particle->Pt(),particle->Eta(),centralityD};
1895
1896          // all genertated primaries including neutral
1897          //if( iMc < stack->GetNprimary() ) {
1898          //fGenTrackMatrix->Fill(vTrackMatrix);
1899          //}
1900
1901          // only charged particles
1902          if(!particle->GetPDG()) continue;
1903          Double_t charge = particle->GetPDG()->Charge()/3.;
1904          if (TMath::Abs(charge) < 0.001)
1905          continue;
1906
1907          // only postive charged
1908          if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.)
1909          continue;
1910
1911          // only negative charged
1912          if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.)
1913          continue;
1914
1915          // physical primary
1916          Bool_t prim = stack->IsPhysicalPrimary(iMc);
1917
1918          // check accepted
1919          if(accCuts->AcceptTrack(particle))
1920          {
1921
1922            if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) {
1923              fGenPrimTrackMatrix->Fill(vTrackMatrix);
1924              Double_t vMCPrimTrackHist[4] = {vtxMC[2],particle->Pt(),particle->Eta(),multMCTrueTracks};
1925              fMCPrimTrackHist->Fill(vMCPrimTrackHist);
1926            }
1927            // fill control histograms
1928            if(fHistogramsOn)
1929              FillHistograms(stack,iMc,AlidNdPtHelper::kAccTracks);
1930
1931            // check multiple found tracks
1932            Int_t multCount = 0;
1933            for(Int_t iRec=0; iRec<multRec; ++iRec)
1934            {
1935              if(iMc == labelsRec[iRec])
1936              {
1937                multCount++;
1938                if(multCount>1)
1939                {
1940                  fRecMultTrackMatrix->Fill(vTrackMatrix);
1941
1942                  // fill control histogram
1943                  if(fHistogramsOn) {
1944                    Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
1945                    Double_t vMCMultRecTrackHist1[3] = {particle->Pt(), particle->Eta(), pid};
1946                    fMCMultRecTrackHist1->Fill(vMCMultRecTrackHist1);
1947                  }
1948                }
1949              }
1950            }
1951
1952            // check reconstructed
1953            for(Int_t iRec=0; iRec<multRec; ++iRec)
1954            {
1955              if(iMc == labelsRec[iRec])
1956              {
1957                fRecTrackMatrix->Fill(vTrackMatrix);
1958
1959                if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) {
1960                  fRecPrimTrackMatrix->Fill(vTrackMatrix);
1961                  //AliESDtrack *track = esdEvent->GetTrack(iRec);
1962                  //if(track && track->GetKinkIndex(0) < 0)
1963                  //  printf("prim kink \n");
1964                  countRecCandle++;
1965                }
1966
1967                if(!prim) fRecSecTrackMatrix->Fill(vTrackMatrix);
1968
1969                // fill control histograms
1970                if(fHistogramsOn)
1971                  FillHistograms(stack,iMc,AlidNdPtHelper::kRecTracks);
1972
1973                break;
1974              }
1975            }
1976          }
1977        }
1978
1979        if(countRecCandle>0) fRecCandleEventMatrix->Fill(vEventMatrix);
1980       }
1981     }
1982   }// end bUseMC
1983
1984   if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;
1985   if(labelsAll) delete [] labelsAll; labelsAll = 0;
1986   if(labelsAcc) delete [] labelsAcc; labelsAcc = 0;
1987   if(labelsRec) delete [] labelsRec; labelsRec = 0;
1988
1989   if(!evtCuts->IsRecVertexRequired() && vtxESD != NULL) delete vtxESD;
1990   //if(trigAna) delete trigAna;
1991
1992 }
1993
1994 //_____________________________________________________________________________
1995 void AlidNdPtAnalysispPb::FillHistograms(TObjArray *const allChargedTracks,Int_t *const labelsAll,Int_t multAll,Int_t *const labelsAcc,Int_t multAcc,Int_t *const labelsRec,Int_t multRec) {
1996  // multiplicity  histograms
1997
1998   if(!allChargedTracks) return;
1999   if(!labelsAll) return;
2000   if(!labelsAcc) return;
2001   if(!labelsRec) return;
2002
2003   Int_t entries = allChargedTracks->GetEntries();
2004   for(Int_t i=0; i<entries; ++i)
2005   {
2006      AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
2007      if(!track) continue;
2008      if(track->Charge() == 0) continue;
2009
2010      Int_t label = TMath::Abs(track->GetLabel());
2011      for(Int_t iAll=0; iAll<multAll; ++iAll) {
2012        if(label == labelsAll[iAll]) {
2013          Double_t v1[2] = {track->Pt(), multAll};
2014          fRecTrackMultHist1[AlidNdPtHelper::kAllTracks]->Fill(v1);
2015        }
2016      }
2017      for(Int_t iAcc=0; iAcc<multAcc; ++iAcc) {
2018        if(label == labelsAcc[iAcc]) {
2019          Double_t v2[2] = {track->Pt(), multAcc};
2020          fRecTrackMultHist1[AlidNdPtHelper::kAccTracks]->Fill(v2);
2021        }
2022      }
2023      for(Int_t iRec=0; iRec<multRec; ++iRec) {
2024        if(label == labelsRec[iRec]) {
2025          Double_t v3[2] = {track->Pt(), multRec};
2026          fRecTrackMultHist1[AlidNdPtHelper::kRecTracks]->Fill(v3);
2027        }
2028      }
2029   }
2030 }
2031
2032 //_____________________________________________________________________________
2033 void AlidNdPtAnalysispPb::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, const Double_t zv, AlidNdPtHelper::TrackObject trackObj, Int_t multRecMult)
2034 {
2035   //
2036   // Fill ESD track and MC histograms
2037   //
2038   if(!esdTrack) return;
2039
2040   Float_t q = esdTrack->Charge();
2041   if(TMath::Abs(q) < 0.001) return;
2042
2043   Float_t pt = esdTrack->Pt();
2044   //Float_t qpt = esdTrack->Pt() * q;
2045   Float_t eta = esdTrack->Eta();
2046   Float_t phi = esdTrack->Phi();
2047
2048   Float_t dca[2], bCov[3];
2049   esdTrack->GetImpactParameters(dca,bCov);
2050
2051   Int_t nClust = esdTrack->GetTPCclusters(0);
2052   Float_t chi2PerCluster = 0.;
2053   if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);
2054
2055
2056   // fill histograms
2057   Double_t values1[3] = {pt,eta,phi};
2058   fRecTrackHist1[trackObj]->Fill(values1);
2059
2060   Double_t values[4] = {zv, pt,eta, multRecMult};
2061   if(trackObj == AlidNdPtHelper::kRecTracks) {
2062     fRecTrackHist->Fill(values);
2063   }
2064
2065   //fill fCentralityTrackHist only if necessary
2066   if (fDimensionsCentralityEstimators > 0) {
2067      fVCentralityTrack[1] = pt;
2068      fVCentralityTrack[2] = eta;
2069      if(trackObj == AlidNdPtHelper::kRecTracks) {
2070        fCentralityTrackHist->Fill(fVCentralityTrack);
2071      }
2072   }
2073
2074   /*
2075   Double_t values2[5] = {nClust,chi2PerCluster,pt,eta,phi};
2076   if(trackObj == AlidNdPtHelper::kRecTracks)
2077   {
2078     if(fHistogramsOn)
2079       fRecTrackHist2->Fill(values2);
2080   }
2081   */
2082
2083   //
2084   // Fill rec vs MC information
2085   //
2086   if(!stack) return;
2087
2088   Int_t label = TMath::Abs(esdTrack->GetLabel());
2089   //if(label == 0) return;
2090
2091   if(label > stack->GetNtrack()) return;
2092   TParticle* particle = stack->Particle(label);
2093   if(!particle) return;
2094
2095   //Bool_t prim = stack->IsPhysicalPrimary(label);
2096   //Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
2097
2098   Int_t motherPdg = -1;
2099   TParticle* mother = 0;
2100
2101   //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
2102   Int_t motherLabel = particle->GetMother(0);
2103   if(motherLabel>0) mother = stack->Particle(motherLabel);
2104   if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
2105   //Int_t mech = particle->GetUniqueID(); // production mechanism
2106
2107   if(!particle->GetPDG()) return;
2108   Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
2109   if(TMath::Abs(gq)<0.001) return;
2110   Float_t gpt = particle->Pt();
2111   Float_t geta = particle->Eta();
2112   //Float_t qgpt = particle->Pt() * gq;
2113   //Float_t gphi = particle->Phi();
2114
2115   Double_t dpt=0;
2116   //printf("pt %f, gpt %f \n",pt,gpt);
2117   if(gpt) dpt = (pt-gpt)/gpt;
2118   Double_t deta = (eta-geta);
2119
2120   // fill histograms
2121   if(trackObj == AlidNdPtHelper::kRecTracks)
2122   {
2123     Double_t vTrackPtCorrelationMatrix[3]={pt,gpt,geta};
2124     fTrackPtCorrelationMatrix->Fill(vTrackPtCorrelationMatrix);
2125
2126     Double_t vRecMCTrackHist1[4]={gpt,geta,dpt,deta};
2127     fRecMCTrackHist1->Fill(vRecMCTrackHist1);
2128   }
2129 }
2130
2131 //_____________________________________________________________________________
2132 void AlidNdPtAnalysispPb::FillHistograms(AliStack *const stack, Int_t label, AlidNdPtHelper::TrackObject trackObj)
2133 {
2134   // Fill MC histograms
2135   if(!stack) return;
2136
2137   if(label > stack->GetNtrack()) return;
2138   TParticle* particle = stack->Particle(label);
2139   if(!particle) return;
2140
2141   Int_t motherPdg = -1;
2142   TParticle* mother = 0;
2143
2144   //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
2145   Int_t motherLabel = particle->GetMother(0);
2146   if(motherLabel>0) mother = stack->Particle(motherLabel);
2147   if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
2148   Int_t mech = particle->GetUniqueID(); // production mechanism
2149
2150   if(!particle->GetPDG()) return;
2151   Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
2152   if(TMath::Abs(gq) < 0.001) return;
2153
2154   Float_t gpt = particle->Pt();
2155   //Float_t qgpt = particle->Pt() * gq;
2156   Float_t geta = particle->Eta();
2157   Float_t gphi = particle->Phi();
2158   //Float_t gpz = particle->Pz();
2159
2160   Bool_t prim = stack->IsPhysicalPrimary(label);
2161   //Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();
2162
2163   Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
2164
2165   //if(prim&&pid==5) printf("pdgcode %d, production mech %d \n",particle->GetPdgCode(),mech);
2166   //if(!prim) printf("motherPdg %d, particle %d, production mech %d\n",motherPdg, particle->GetPdgCode(),mech);
2167
2168   //
2169   // fill histogram
2170   //
2171   Double_t vMCTrackHist1[3] = {gpt,geta,gphi};
2172   fMCTrackHist1[trackObj]->Fill(vMCTrackHist1);
2173
2174   Double_t vMCPrimTrackHist1[5] = {gpt,geta,pid,mech,motherPdg};
2175   Double_t vMCPrimTrackHist2[5] = {TMath::Abs(particle->GetPdgCode()),mech,motherPdg};
2176   //if(prim && AliPWG0Helper::IsPrimaryCharged(particle, label)) fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
2177   if(prim) {
2178     fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
2179     if(pid == 5) fMCPrimTrackHist2[trackObj]->Fill(vMCPrimTrackHist2);
2180   }
2181   else {
2182     fMCSecTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
2183   }
2184
2185 }
2186
2187 //_____________________________________________________________________________
2188 Long64_t AlidNdPtAnalysispPb::Merge(TCollection* const list)
2189 {
2190   //  init if not done already
2191   if (!fIsInit) { Init(); }
2192
2193   // Merge list of objects (needed by PROOF)
2194
2195   if (!list)
2196   return 0;
2197
2198   if (list->IsEmpty())
2199   return 1;
2200
2201   TIterator* iter = list->MakeIterator();
2202   TObject* obj = 0;
2203
2204   //
2205   //TList *collPhysSelection = new TList;
2206
2207   // collection of generated histograms
2208
2209   Int_t count=0;
2210   while((obj = iter->Next()) != 0) {
2211     AlidNdPtAnalysispPb* entry = dynamic_cast<AlidNdPtAnalysispPb*>(obj);
2212     if (entry == 0) continue;
2213
2214     // physics selection
2215     //printf("entry->GetPhysicsTriggerSelection() %p \n", entry->GetPhysicsTriggerSelection());
2216     //collPhysSelection->Add(entry->GetPhysicsTriggerSelection());
2217
2218     //
2219     fRecEventHist->Add(entry->fRecEventHist);
2220     fRecTrackHist->Add(entry->fRecTrackHist);
2221     fEventCount->Add(entry->fEventCount);
2222     fEventMultHist->Add(entry->fEventMultHist);
2223     if (fDimensionsCentralityEstimators > 0) {
2224         fCentralityEventHist->Add(entry->fCentralityEventHist);
2225         fCentralityTrackHist->Add(entry->fCentralityTrackHist);
2226     }
2227
2228     //
2229     fEventMultCorrelationMatrix->Add(entry->fEventMultCorrelationMatrix);
2230     fTrackPtCorrelationMatrix->Add(entry->fTrackPtCorrelationMatrix);
2231
2232     //
2233     fGenEventMatrix->Add(entry->fGenEventMatrix);
2234     fGenSDEventMatrix->Add(entry->fGenSDEventMatrix);
2235     fGenDDEventMatrix->Add(entry->fGenDDEventMatrix);
2236     fGenNDEventMatrix->Add(entry->fGenNDEventMatrix);
2237     fGenNSDEventMatrix->Add(entry->fGenNSDEventMatrix);
2238
2239     fTriggerEventMatrix->Add(entry->fTriggerEventMatrix);
2240     fTriggerSDEventMatrix->Add(entry->fTriggerSDEventMatrix);
2241     fTriggerDDEventMatrix->Add(entry->fTriggerDDEventMatrix);
2242     fTriggerNDEventMatrix->Add(entry->fTriggerNDEventMatrix);
2243     fTriggerNSDEventMatrix->Add(entry->fTriggerNSDEventMatrix);
2244
2245     fRecEventMatrix->Add(entry->fRecEventMatrix);
2246     fRecSDEventMatrix->Add(entry->fRecSDEventMatrix);
2247     fRecDDEventMatrix->Add(entry->fRecDDEventMatrix);
2248     fRecNDEventMatrix->Add(entry->fRecNDEventMatrix);
2249     fRecNSDEventMatrix->Add(entry->fRecNSDEventMatrix);
2250
2251     fRecCandleEventMatrix->Add(entry->fRecCandleEventMatrix);
2252     //
2253     fGenTrackEventMatrix->Add(entry->fGenTrackEventMatrix);
2254     fGenTrackSDEventMatrix->Add(entry->fGenTrackSDEventMatrix);
2255     fGenTrackDDEventMatrix->Add(entry->fGenTrackDDEventMatrix);
2256     fGenTrackNDEventMatrix->Add(entry->fGenTrackNDEventMatrix);
2257     fGenTrackNSDEventMatrix->Add(entry->fGenTrackNSDEventMatrix);
2258
2259     fTriggerTrackEventMatrix->Add(entry->fTriggerTrackEventMatrix);
2260     fTriggerTrackSDEventMatrix->Add(entry->fTriggerTrackSDEventMatrix);
2261     fTriggerTrackDDEventMatrix->Add(entry->fTriggerTrackDDEventMatrix);
2262     fTriggerTrackNDEventMatrix->Add(entry->fTriggerTrackNDEventMatrix);
2263     fTriggerTrackNSDEventMatrix->Add(entry->fTriggerTrackNSDEventMatrix);
2264
2265     fRecTrackEventMatrix->Add(entry->fRecTrackEventMatrix);
2266     fRecTrackSDEventMatrix->Add(entry->fRecTrackSDEventMatrix);
2267     fRecTrackDDEventMatrix->Add(entry->fRecTrackDDEventMatrix);
2268     fRecTrackNDEventMatrix->Add(entry->fRecTrackNDEventMatrix);
2269     fRecTrackNSDEventMatrix->Add(entry->fRecTrackNSDEventMatrix);
2270
2271     //
2272     fGenTrackMatrix->Add(entry->fGenTrackMatrix);
2273     fGenPrimTrackMatrix->Add(entry->fGenPrimTrackMatrix);
2274     fRecPrimTrackMatrix->Add(entry->fRecPrimTrackMatrix);
2275     //
2276     fRecTrackMatrix->Add(entry->fRecTrackMatrix);
2277     fRecSecTrackMatrix->Add(entry->fRecSecTrackMatrix);
2278     //
2279     fRecMultTrackMatrix->Add(entry->fRecMultTrackMatrix);
2280     fMCPrimTrackHist->Add(entry->fMCPrimTrackHist);
2281
2282     //
2283     // control analysis histograms
2284     //
2285     fMCEventHist1->Add(entry->fMCEventHist1);
2286     fRecEventHist1->Add(entry->fRecEventHist1);
2287     fRecEventHist2->Add(entry->fRecEventHist2);
2288     fRecMCEventHist1->Add(entry->fRecMCEventHist1);
2289     fRecMCEventHist2->Add(entry->fRecMCEventHist2);
2290     fRecMCEventHist3->Add(entry->fRecMCEventHist3);
2291
2292     for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
2293       fMCTrackHist1[i]->Add(entry->fMCTrackHist1[i]);
2294
2295       fMCPrimTrackHist1[i]->Add(entry->fMCPrimTrackHist1[i]);
2296       fMCPrimTrackHist2[i]->Add(entry->fMCPrimTrackHist2[i]);
2297       fMCSecTrackHist1[i]->Add(entry->fMCSecTrackHist1[i]);
2298
2299       fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);
2300       fRecTrackMultHist1[i]->Add(entry->fRecTrackMultHist1[i]);
2301     }
2302     fRecMCTrackHist1->Add(entry->fRecMCTrackHist1);
2303     fMCMultRecTrackHist1->Add(entry->fMCMultRecTrackHist1);
2304     fRecTrackHist2->Add(entry->fRecTrackHist2);
2305
2306
2307   count++;
2308   }
2309
2310   //AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();
2311   //trigSelection->Merge(collPhysSelection);
2312   //if(collPhysSelection) delete collPhysSelection;
2313
2314 return count;
2315 }
2316
2317 //____________________________________________________________________
2318 Bool_t AlidNdPtAnalysispPb::IsBinZeroTrackSPDvtx(const AliESDEvent* esdEvent) {
2319 //
2320 // check 0-bin
2321 // for LHC background calculation
2322 // return kTRUE if vertex not reconstructed or
2323 // track multiplicity == 0
2324 //
2325 if(!esdEvent) return kFALSE;
2326
2327   // check vertex
2328   const AliESDVertex *vertex = esdEvent->GetPrimaryVertexTracks();
2329   if(!vertex) return kTRUE;
2330
2331   if(vertex->GetNContributors() < 1) {
2332     // SPD vertex
2333     vertex = esdEvent->GetPrimaryVertexSPD();
2334     if(!vertex) return kTRUE;
2335   }
2336   if(TMath::Abs(vertex->GetZv()) > 15.0) return kTRUE;
2337   if( !AlidNdPtHelper::TestRecVertex(vertex, esdEvent->GetPrimaryVertexSPD(), AlidNdPtHelper::kTPCITSHybridTrackSPDvtx) )
2338     return kTRUE;
2339
2340 return kFALSE;
2341 }
2342
2343 //____________________________________________________________________
2344 Bool_t AlidNdPtAnalysispPb::IsBinZeroSPDvtx(const AliESDEvent* esdEvent) {
2345 //
2346 // check 0-bin
2347 // for LHC background calculation
2348 // return kTRUE if vertex not reconstructed or
2349 // tracklet multiplicity == 0
2350 //
2351 if(!esdEvent) return kFALSE;
2352
2353   // check vertex
2354   const AliESDVertex* vertex = esdEvent->GetPrimaryVertexSPD();
2355   if(!vertex) return kTRUE;
2356   if(TMath::Abs(vertex->GetZv()) > 15.0) return kTRUE;
2357   if( !AlidNdPtHelper::TestRecVertex(vertex, esdEvent->GetPrimaryVertexSPD(), AlidNdPtHelper::kTPCSPDvtx) ) return kTRUE;
2358
2359 return kFALSE;
2360 }
2361
2362 //_____________________________________________________________________________
2363 void AlidNdPtAnalysispPb::Analyse()
2364 {
2365   //  init if not done already
2366   if (!fIsInit) { Init(); }
2367
2368   // Analyse histograms
2369   //
2370   TH1::AddDirectory(kFALSE);
2371   TH1 *h=0, *h1=0, *h2=0, *h2c = 0;
2372   THnSparse *hs=0;
2373   TH2 *h2D=0;
2374
2375   char name[256];
2376   TObjArray *aFolderObj = new TObjArray;
2377   if(!aFolderObj) return;
2378
2379   //
2380   // LHC backgraund in all and 0-bins
2381   //
2382   //AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection();
2383   //trigSel->SaveHistograms("physics_selection");
2384
2385   //
2386   // Reconstructed event vertex
2387   //
2388   h = fRecEventHist1->Projection(0);
2389   h->SetName("Xv");
2390   aFolderObj->Add(h);
2391
2392   h = fRecEventHist1->Projection(1);
2393   h->SetName("Yv");
2394   aFolderObj->Add(h);
2395
2396   h = fRecEventHist1->Projection(2);
2397   h->SetName("Zv");
2398   aFolderObj->Add(h);
2399
2400   //
2401   // multiplicity
2402   //
2403   h = fRecEventHist2->Projection(1);
2404   h->SetName("multMB");
2405   aFolderObj->Add(h);
2406
2407   h = fRecEventHist2->Projection(2);
2408   h->SetName("multiplicity");
2409   aFolderObj->Add(h);
2410
2411   h2D = fRecEventHist2->Projection(0,1);
2412   h2D->SetName("Zv_vs_multiplicity_MB");
2413   aFolderObj->Add(h2D);
2414
2415   //
2416   // reconstructed pt histograms
2417   //
2418   h = fRecTrackHist1[0]->Projection(0);
2419   h->Scale(1.,"width");
2420   h->SetName("pt_all_ch");
2421   aFolderObj->Add(h);
2422
2423   h = fRecTrackHist1[1]->Projection(0);
2424   h->Scale(1.,"width");
2425   h->SetName("pt_acc");
2426   aFolderObj->Add(h);
2427
2428   h = fRecTrackHist1[2]->Projection(0);
2429   h->Scale(1.,"width");
2430   h->SetName("pt_rec");
2431   aFolderObj->Add(h);
2432
2433   //
2434   // reconstructed eta histograms
2435   //
2436   h = fRecTrackHist1[0]->Projection(1);
2437   h->SetName("eta_all_ch");
2438   aFolderObj->Add(h);
2439
2440   h = fRecTrackHist1[1]->Projection(1);
2441   h->SetName("eta_acc");
2442   aFolderObj->Add(h);
2443
2444   h = fRecTrackHist1[2]->Projection(1);
2445   h->SetName("eta_rec");
2446   aFolderObj->Add(h);
2447
2448   //
2449   // reconstructed phi histograms
2450   //
2451   h = fRecTrackHist1[0]->Projection(2);
2452   h->SetName("phi_all_ch");
2453   aFolderObj->Add(h);
2454
2455   h = fRecTrackHist1[1]->Projection(2);
2456   h->SetName("phi_acc");
2457   aFolderObj->Add(h);
2458
2459   h = fRecTrackHist1[2]->Projection(2);
2460   h->SetName("phi_rec");
2461   aFolderObj->Add(h);
2462
2463   //
2464   // reconstructed eta:pt histograms
2465   //
2466   h2D = fRecTrackHist1[0]->Projection(1,0);
2467   h2D->SetName("pt_eta_all_ch");
2468   aFolderObj->Add(h2D);
2469
2470   h2D = fRecTrackHist1[1]->Projection(1,0);
2471   h2D->SetName("pt_eta_acc");
2472   aFolderObj->Add(h2D);
2473
2474   h2D = fRecTrackHist1[2]->Projection(1,0);
2475   h2D->SetName("pt_eta_rec");
2476   aFolderObj->Add(h2D);
2477
2478   //
2479   // reconstructed phi:pt histograms
2480   //
2481   h2D = fRecTrackHist1[0]->Projection(2,0);
2482   h2D->SetName("pt_phi_all_ch");
2483   aFolderObj->Add(h2D);
2484
2485   h2D = fRecTrackHist1[1]->Projection(2,0);
2486   h2D->SetName("pt_phi_acc");
2487   aFolderObj->Add(h2D);
2488
2489   h2D = fRecTrackHist1[2]->Projection(2,0);
2490   h2D->SetName("pt_phi_rec");
2491   aFolderObj->Add(h2D);
2492
2493   //
2494   // reconstructed phi:eta histograms
2495   //
2496   h2D = fRecTrackHist1[0]->Projection(2,1);
2497   h2D->SetName("eta_phi_all_ch");
2498   aFolderObj->Add(h2D);
2499
2500   h2D = fRecTrackHist1[1]->Projection(2,1);
2501   h2D->SetName("eta_phi_acc");
2502   aFolderObj->Add(h2D);
2503
2504   h2D = fRecTrackHist1[2]->Projection(2,1);
2505   h2D->SetName("eta_phi_rec");
2506   aFolderObj->Add(h2D);
2507
2508   //
2509   // reconstructed nClust, chi2 vs pt, eta, phi
2510   //
2511   if(fHistogramsOn) {
2512
2513     h2D = fRecTrackHist2->Projection(0,1);
2514     h2D->SetName("nClust_chi2_rec");
2515     aFolderObj->Add(h2D);
2516
2517     h2D = fRecTrackHist2->Projection(0,2);
2518     h2D->SetName("nClust_pt_rec");
2519     aFolderObj->Add(h2D);
2520
2521     h2D = fRecTrackHist2->Projection(0,3);
2522     h2D->SetName("nClust_eta_rec");
2523     aFolderObj->Add(h2D);
2524
2525     h2D = fRecTrackHist2->Projection(0,4);
2526     h2D->SetName("nClust_phi_rec");
2527     aFolderObj->Add(h2D);
2528
2529     h2D = fRecTrackHist2->Projection(1,2);
2530     h2D->SetName("chi2_pt_rec");
2531     aFolderObj->Add(h2D);
2532
2533     h2D = fRecTrackHist2->Projection(1,3);
2534     h2D->SetName("chi2_eta_rec");
2535     aFolderObj->Add(h2D);
2536
2537     h2D = fRecTrackHist2->Projection(1,4);
2538     h2D->SetName("chi2_phi_rec");
2539     aFolderObj->Add(h2D);
2540
2541   }
2542
2543   //
2544   // calculate corrections for empty events
2545   // with multMB==0
2546   //
2547
2548   //
2549   // normalised zv to generate zv for triggered events
2550   //
2551   h = fRecEventHist2->Projection(0);
2552   if( h->Integral() ) h->Scale(1./h->Integral());
2553   h->SetName("zv_distribution_norm");
2554   aFolderObj->Add(h);
2555
2556   //
2557   // MC available
2558   //
2559   if(IsUseMCInfo()) {
2560
2561   //
2562   // Event vertex resolution
2563   //
2564   h2D = fRecMCEventHist2->Projection(0,2);
2565   h2D->SetName("DeltaXv_vs_mult");
2566   aFolderObj->Add(h2D);
2567
2568   h2D = fRecMCEventHist2->Projection(1,2);
2569   h2D->SetName("DeltaZv_vs_mult");
2570   aFolderObj->Add(h2D);
2571
2572   //
2573   // normalised zv to get trigger/trigger+vertex event differences
2574   // F(zv) = E_trig(zv,0)/Int(E_trig(zv,0) / Sum(E_trigvtx(zv,n))/Sum(Int(E_trigvtx(zv,n))dzv)
2575   //
2576   fTriggerEventMatrix->GetAxis(1)->SetRangeUser(0.,0.);
2577   h = fTriggerEventMatrix->Projection(0);
2578   h2D = fTriggerEventMatrix->Projection(0,1);
2579   if(h2D->Integral()) h->Scale(1./h2D->Integral());
2580
2581   h1 = fRecEventMatrix->Projection(0);
2582   h2D = fRecEventMatrix->Projection(0,1);
2583   if(h2D->Integral()) h1->Scale(1./h2D->Integral());
2584
2585   h->Divide(h1);
2586   h->SetName("zv_empty_events_norm");
2587   aFolderObj->Add(h);
2588
2589   fTriggerEventMatrix->GetAxis(1)->SetRange(1,fTriggerEventMatrix->GetAxis(1)->GetNbins());
2590
2591   //
2592   // rec. vs true multiplicity correlation matrix
2593   //
2594   hs = (THnSparse*)fEventMultCorrelationMatrix->Clone("event_mult_correlation_matrix");
2595   aFolderObj->Add(hs);
2596
2597   //
2598   // rec. vs true track pt correlation matrix
2599   //
2600   hs = (THnSparse*)fTrackPtCorrelationMatrix->Clone("track_pt_correlation_matrix");
2601   aFolderObj->Add(hs);
2602
2603   //
2604   // trigger efficiency for INEL
2605   //
2606   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fGenEventMatrix->Projection(0),"zv_trig_INEL_eff_matrix");
2607   aFolderObj->Add(h);
2608
2609   //
2610   // trigger efficiency for NSD
2611   //
2612   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerNSDEventMatrix->Projection(0),fGenNSDEventMatrix->Projection(0),"zv_trig_NSD_eff_matrix");
2613   aFolderObj->Add(h);
2614
2615   //
2616   // trigger bias correction (MB to ND)
2617   //
2618   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenNDEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoND_corr_matrix");
2619   aFolderObj->Add(hs);
2620
2621   h = AlidNdPtHelper::GenerateCorrMatrix(fGenNDEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoND_corr_matrix");
2622   aFolderObj->Add(h);
2623
2624
2625   h = AlidNdPtHelper::GenerateCorrMatrix(fGenNDEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoND_corr_matrix");
2626
2627   aFolderObj->Add(h);
2628
2629   //
2630   // trigger bias correction (MB to NSD)
2631   //
2632   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenNSDEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoNSD_corr_matrix");
2633   aFolderObj->Add(hs);
2634
2635   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenNSDEventMatrix->Projection(0,1),fTriggerEventMatrix->Projection(0,1),"zv_mult_trig_MBtoNSD_corr_matrix_2D");
2636   aFolderObj->Add(h2D);
2637
2638   h = AlidNdPtHelper::GenerateCorrMatrix(fGenNSDEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoNSD_corr_matrix");
2639   aFolderObj->Add(h);
2640
2641
2642   h = AlidNdPtHelper::GenerateCorrMatrix(fGenNSDEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoNSD_corr_matrix");
2643   aFolderObj->Add(h);
2644
2645
2646   //
2647   // trigger bias correction (MB to INEL)
2648   //
2649   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoInel_corr_matrix");
2650   aFolderObj->Add(hs);
2651
2652   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0,1),fTriggerEventMatrix->Projection(0,1),"zv_mult_trig_MBtoInel_corr_matrix_2D");
2653   aFolderObj->Add(h2D);
2654
2655   h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoInel_corr_matrix");
2656   aFolderObj->Add(h);
2657
2658   h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoInel_corr_matrix");
2659   aFolderObj->Add(h);
2660
2661
2662   //
2663   // event vertex reconstruction correction (MB)
2664   //
2665   hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix,fRecEventMatrix,"zv_mult_event_corr_matrix");
2666   aFolderObj->Add(hs);
2667
2668   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0,1),fRecEventMatrix->Projection(0,1),"zv_mult_event_corr_matrix_2D");
2669   aFolderObj->Add(h2D);
2670
2671
2672   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(1),fRecEventMatrix->Projection(1),"mult_event_corr_matrix");
2673   aFolderObj->Add(h);
2674
2675
2676   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fRecEventMatrix->Projection(0),"zv_event_corr_matrix");
2677   aFolderObj->Add(h);
2678
2679   //
2680   // track-event trigger bias correction (MB to ND)
2681   //
2682
2683   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNDEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoND_corr_matrix");
2684   aFolderObj->Add(hs);
2685
2686   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNDEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoND_corr_matrix");
2687   aFolderObj->Add(h2D);
2688
2689   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNDEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoND_corr_matrix");
2690   aFolderObj->Add(h2D);
2691
2692   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNDEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoND_corr_matrix");
2693   aFolderObj->Add(h2D);
2694
2695   //
2696   // track-event trigger bias correction (MB to NSD)
2697   //
2698   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNSDEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoNSD_corr_matrix");
2699   aFolderObj->Add(hs);
2700
2701   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNSDEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoNSD_corr_matrix");
2702   aFolderObj->Add(h2D);
2703
2704   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNSDEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoNSD_corr_matrix");
2705   aFolderObj->Add(h2D);
2706
2707   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackNSDEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoNSD_corr_matrix");
2708   aFolderObj->Add(h2D);
2709
2710
2711   //
2712   // track-event trigger bias correction (MB to INEL)
2713   //
2714   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoInel_corr_matrix");
2715   aFolderObj->Add(hs);
2716
2717   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoInel_corr_matrix");
2718   aFolderObj->Add(h2D);
2719
2720   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoInel_corr_matrix");
2721   aFolderObj->Add(h2D);
2722
2723   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoInel_corr_matrix");
2724   aFolderObj->Add(h2D);
2725
2726   // efficiency
2727
2728   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1),fGenTrackEventMatrix->Projection(1),"pt_track_trig_MBtoInel_eff_matrix");
2729   aFolderObj->Add(h);
2730
2731
2732   //
2733   // track-event vertex reconstruction correction (MB)
2734   //
2735   hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix,fRecTrackEventMatrix,"zv_pt_eta_track_event_corr_matrix");
2736   aFolderObj->Add(hs);
2737
2738   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,2),fRecTrackEventMatrix->Projection(1,2),"eta_pt_track_event_corr_matrix");
2739   aFolderObj->Add(h2D);
2740
2741   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,0),fRecTrackEventMatrix->Projection(1,0),"pt_zv_track_event_corr_matrix");
2742   aFolderObj->Add(h2D);
2743
2744   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(2,0),fRecTrackEventMatrix->Projection(2,0),"zv_eta_track_event_corr_matrix");
2745   aFolderObj->Add(h2D);
2746
2747   // efficiency
2748
2749   h = AlidNdPtHelper::GenerateCorrMatrix(fRecTrackEventMatrix->Projection(1),fTriggerTrackEventMatrix->Projection(1),"pt_track_event_eff_matrix");
2750   aFolderObj->Add(h);
2751
2752
2753   //
2754   // track rec. efficiency correction
2755   //
2756   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix,fRecPrimTrackMatrix,"zv_pt_eta_track_corr_matrix");
2757   aFolderObj->Add(hs);
2758
2759   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,2),fRecPrimTrackMatrix->Projection(1,2),"eta_pt_track_corr_matrix");
2760   aFolderObj->Add(h2D);
2761
2762   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,0),fRecPrimTrackMatrix->Projection(1,0),"pt_zv_track_corr_matrix");
2763   aFolderObj->Add(h2D);
2764
2765   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2,0),fRecPrimTrackMatrix->Projection(2,0),"zv_eta_track_corr_matrix");
2766   aFolderObj->Add(h2D);
2767
2768
2769   h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(0),fRecPrimTrackMatrix->Projection(0),"zv_track_corr_matrix");
2770   aFolderObj->Add(h);
2771
2772   h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1),fRecPrimTrackMatrix->Projection(1),"pt_track_corr_matrix");
2773   aFolderObj->Add(h);
2774
2775   // efficiency
2776
2777   h = AlidNdPtHelper::GenerateCorrMatrix(fRecPrimTrackMatrix->Projection(1), fGenPrimTrackMatrix->Projection(1),"pt_track_eff_matrix");
2778   aFolderObj->Add(h);
2779
2780   h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2),fRecPrimTrackMatrix->Projection(2),"eta_track_corr_matrix");
2781   aFolderObj->Add(h);
2782
2783   //
2784   // secondary track contamination correction
2785   //
2786   //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");
2787   hs = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");
2788   aFolderObj->Add(hs);
2789
2790   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_track_cont_matrix");
2791   aFolderObj->Add(h2D);
2792
2793   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_track_cont_matrix");
2794   aFolderObj->Add(h2D);
2795
2796   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_track_cont_matrix");
2797   aFolderObj->Add(h2D);
2798
2799   h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_track_cont_matrix");
2800   aFolderObj->Add(h);
2801
2802
2803   h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_track_cont_matrix");
2804   aFolderObj->Add(h);
2805
2806   h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_track_cont_matrix");
2807   aFolderObj->Add(h);
2808
2809   //
2810   // multiple track reconstruction correction
2811   //
2812   //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");
2813   hs = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");
2814   aFolderObj->Add(hs);
2815
2816   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_mult_track_cont_matrix");
2817   aFolderObj->Add(h2D);
2818
2819   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_mult_track_cont_matrix");
2820   aFolderObj->Add(h2D);
2821
2822   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_mult_track_cont_matrix");
2823   aFolderObj->Add(h2D);
2824
2825   h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_mult_track_cont_matrix");
2826   aFolderObj->Add(h);
2827
2828   h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_mult_track_cont_matrix");
2829   aFolderObj->Add(h);
2830
2831   h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_mult_track_cont_matrix");
2832   aFolderObj->Add(h);
2833
2834   //
2835   // Control histograms
2836   //
2837
2838   if(fHistogramsOn) {
2839
2840   // Efficiency electrons, muons, pions, kaons, protons, all
2841   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,1);
2842   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,1);
2843   h1 = fMCPrimTrackHist1[1]->Projection(0);
2844   h2 = fMCPrimTrackHist1[2]->Projection(0);
2845   h2c = (TH1D *)h2->Clone();
2846   h2c->Divide(h1);
2847   h2c->SetName("eff_pt_electrons");
2848   aFolderObj->Add(h2c);
2849
2850   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(2,2);
2851   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(2,2);
2852   h1 = fMCPrimTrackHist1[1]->Projection(0);
2853   h2 = fMCPrimTrackHist1[2]->Projection(0);
2854   h2c = (TH1D *)h2->Clone();
2855   h2c->Divide(h1);
2856   h2c->SetName("eff_pt_muons");
2857   aFolderObj->Add(h2c);
2858
2859   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(3,3);
2860   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(3,3);
2861   h1 = fMCPrimTrackHist1[1]->Projection(0);
2862   h2 = fMCPrimTrackHist1[2]->Projection(0);
2863   h2c = (TH1D *)h2->Clone();
2864   h2c->Divide(h1);
2865   h2c->SetName("eff_pt_pions");
2866   aFolderObj->Add(h2c);
2867
2868   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(4,4);
2869   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(4,4);
2870   h1 = fMCPrimTrackHist1[1]->Projection(0);
2871   h2 = fMCPrimTrackHist1[2]->Projection(0);
2872   h2c = (TH1D *)h2->Clone();
2873   h2c->Divide(h1);
2874   h2c->SetName("eff_pt_kaons");
2875   aFolderObj->Add(h2c);
2876
2877   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(5,5);
2878   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(5,5);
2879   h1 = fMCPrimTrackHist1[1]->Projection(0);
2880   h2 = fMCPrimTrackHist1[2]->Projection(0);
2881   h2c = (TH1D *)h2->Clone();
2882   h2c->Divide(h1);
2883   h2c->SetName("eff_pt_protons");
2884   aFolderObj->Add(h2c);
2885
2886   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,5);
2887   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,5);
2888   h1 = fMCPrimTrackHist1[1]->Projection(0);
2889   h2 = fMCPrimTrackHist1[2]->Projection(0);
2890   h2c = (TH1D *)h2->Clone();
2891   h2c->Divide(h1);
2892   h2c->SetName("eff_pt_selected");
2893   aFolderObj->Add(h2c);
2894
2895   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,6);
2896   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,6);
2897   h1 = fMCPrimTrackHist1[1]->Projection(0);
2898   h2 = fMCPrimTrackHist1[2]->Projection(0);
2899   h2c = (TH1D *)h2->Clone();
2900   h2c->Divide(h1);
2901   h2c->SetName("eff_pt_all");
2902   aFolderObj->Add(h2c);
2903
2904   fMCPrimTrackHist1[1]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[1]->GetAxis(1)->GetNbins());
2905   fMCPrimTrackHist1[2]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[2]->GetAxis(1)->GetNbins());
2906
2907   //  pt spetra
2908   // - rec, primaries, secondaries
2909   // - primaries (pid)
2910   // - secondaries (pid)
2911   // - secondaries (mech)
2912   // - secondaries (mother)
2913   //
2914
2915   TH1D *mcPtAccall = fMCTrackHist1[1]->Projection(0);
2916   mcPtAccall->SetName("mc_pt_acc_all");
2917   aFolderObj->Add(mcPtAccall);
2918
2919   TH1D *mcPtAccprim = fMCPrimTrackHist1[1]->Projection(0);
2920   mcPtAccprim->SetName("mc_pt_acc_prim");
2921   aFolderObj->Add(mcPtAccprim);
2922
2923   TH1D *mcPtRecall = fMCTrackHist1[2]->Projection(0);
2924   mcPtRecall->SetName("mc_pt_rec_all");
2925   aFolderObj->Add(mcPtRecall);
2926
2927   TH1D *mcPtRecprim = fMCPrimTrackHist1[2]->Projection(0);
2928   mcPtRecprim->SetName("mc_pt_rec_prim");
2929   aFolderObj->Add(mcPtRecprim);
2930
2931   TH1D *mcPtRecsec = fMCSecTrackHist1[2]->Projection(0);
2932   mcPtRecsec->SetName("mc_pt_rec_sec");
2933   aFolderObj->Add(mcPtRecsec);
2934
2935   for(Int_t i = 0; i<6; i++)
2936   {
2937     snprintf(name,256,"mc_pt_rec_prim_pid_%d",i);
2938     fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
2939     h = fMCPrimTrackHist1[2]->Projection(0);
2940     h->SetName(name);
2941     aFolderObj->Add(h);
2942
2943     snprintf(name,256,"mc_pt_rec_sec_pid_%d",i);
2944     fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
2945     h = fMCSecTrackHist1[2]->Projection(0);
2946     h->SetName(name);
2947     aFolderObj->Add(h);
2948
2949     // production mechanisms for given pid
2950     fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
2951
2952     for(Int_t j=0; j<20; j++) {
2953       if(j == 4) {
2954         // decay
2955
2956         snprintf(name,256,"mc_pt_rec_sec_pid_%d_decay",i);
2957         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
2958         h = fMCSecTrackHist1[2]->Projection(0);
2959         h->SetName(name);
2960         aFolderObj->Add(h);
2961
2962         snprintf(name,256,"mc_eta_rec_sec_pid_%d_decay",i);
2963         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
2964         h = fMCSecTrackHist1[2]->Projection(1);
2965         h->SetName(name);
2966         aFolderObj->Add(h);
2967
2968         snprintf(name,256,"mc_mother_rec_sec_pid_%d_decay",i);
2969         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
2970         h = fMCSecTrackHist1[2]->Projection(4);
2971         h->SetName(name);
2972         aFolderObj->Add(h);
2973
2974       } else if (j == 5) {
2975        // conversion
2976
2977         snprintf(name,256,"mc_pt_rec_sec_pid_%d_conv",i);
2978         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
2979         h = fMCSecTrackHist1[2]->Projection(0);
2980         h->SetName(name);
2981         aFolderObj->Add(h);
2982
2983         snprintf(name,256,"mc_eta_rec_sec_pid_%d_conv",i);
2984         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
2985         h = fMCSecTrackHist1[2]->Projection(1);
2986         h->SetName(name);
2987         aFolderObj->Add(h);
2988
2989         snprintf(name,256,"mc_mother_rec_sec_pid_%d_conv",i);
2990         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
2991         h = fMCSecTrackHist1[2]->Projection(4);
2992         h->SetName(name);
2993         aFolderObj->Add(h);
2994
2995      } else if (j == 13) {
2996        // mat
2997
2998         snprintf(name,256,"mc_pt_rec_sec_pid_%d_mat",i);
2999         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
3000         h = fMCSecTrackHist1[2]->Projection(0);
3001         h->SetName(name);
3002         aFolderObj->Add(h);
3003
3004         snprintf(name,256,"mc_eta_rec_sec_pid_%d_mat",i);
3005         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
3006         h = fMCSecTrackHist1[2]->Projection(1);
3007         h->SetName(name);
3008         aFolderObj->Add(h);
3009
3010         snprintf(name,256,"mc_eta_mother_rec_sec_pid_%d_mat",i);
3011         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
3012         h = fMCSecTrackHist1[2]->Projection(4,1);
3013         h->SetName(name);
3014         aFolderObj->Add(h);
3015
3016         snprintf(name,256,"mc_mother_rec_sec_pid_%d_mat",i);
3017         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
3018         h = fMCSecTrackHist1[2]->Projection(4);
3019         h->SetName(name);
3020         aFolderObj->Add(h);
3021
3022         snprintf(name,256,"mc_pt_mother_rec_sec_pid_%d_mat",i);
3023         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
3024         h = fMCSecTrackHist1[2]->Projection(4,0);
3025         h->SetName(name);
3026         aFolderObj->Add(h);
3027
3028      } else {
3029        continue;
3030      }
3031    }
3032
3033   }
3034   } // end fHistogramOn
3035
3036   //
3037   //  resolution histograms
3038   //  only for reconstructed tracks
3039   //
3040
3041   TH2F *h2F=0;
3042   TCanvas * c = new TCanvas("resol","resol");
3043   c->cd();
3044
3045   //
3046   fRecMCTrackHist1->GetAxis(1)->SetRangeUser(-0.8,0.79);
3047
3048   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);
3049   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
3050   h->SetXTitle("p_{tmc} (GeV/c)");
3051   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");
3052   h->Draw();
3053   h->SetName("pt_resolution_vs_mcpt");
3054   aFolderObj->Add(h);
3055
3056   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);
3057   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
3058   h->SetXTitle("p_{tmc} (GeV/c)");
3059   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");
3060   h->Draw();
3061   h->SetName("dpt_mean_vs_mcpt");
3062   aFolderObj->Add(h);
3063
3064   //
3065   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);
3066   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
3067   h->SetXTitle("p_{tmc} (GeV/c)");
3068   h->SetYTitle("(#eta-#eta_{mc}) resolution");
3069   h->Draw();
3070   h->SetName("eta_resolution_vs_mcpt");
3071   aFolderObj->Add(h);
3072
3073   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);
3074   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
3075   h->SetXTitle("p_{tmc} (GeV/c)");
3076   h->SetYTitle("(#eta-mc#eta) mean");
3077   h->Draw();
3078   h->SetName("deta_mean_vs_mcpt");
3079   aFolderObj->Add(h);
3080
3081   //
3082   fRecMCTrackHist1->GetAxis(1)->SetRange(1,fRecMCTrackHist1->GetAxis(1)->GetNbins());
3083
3084   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);
3085   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
3086   h->SetXTitle("#eta_{mc}");
3087   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");
3088   h->Draw();
3089   h->SetName("pt_resolution_vs_mceta");
3090   aFolderObj->Add(h);
3091
3092   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);
3093   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
3094   h->SetXTitle("#eta_{mc}");
3095   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");
3096   h->Draw();
3097   h->SetName("dpt_mean_vs_mceta");
3098   aFolderObj->Add(h);
3099
3100   //
3101   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);
3102   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
3103   h->SetXTitle("#eta_{mc}");
3104   h->SetYTitle("(#eta-#eta_{mc}) resolution");
3105   h->Draw();
3106   h->SetName("eta_resolution_vs_mceta");
3107   aFolderObj->Add(h);
3108
3109   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);
3110   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
3111   h->SetXTitle("#eta_{mc}");
3112   h->SetYTitle("(#eta-mc#eta) mean");
3113   h->Draw();
3114   h->SetName("deta_mean_vs_mceta");
3115   aFolderObj->Add(h);
3116
3117   fRecMCTrackHist1->GetAxis(0)->SetRange(1,fRecMCTrackHist1->GetAxis(0)->GetNbins());
3118
3119   } // end use MC info
3120
3121   // export objects to analysis folder
3122   fAnalysisFolder = ExportToFolder(aFolderObj);
3123
3124   // delete only TObjArray
3125   if(aFolderObj) delete aFolderObj;
3126 }
3127
3128 //_____________________________________________________________________________
3129 TFolder* AlidNdPtAnalysispPb::ExportToFolder(TObjArray * const array)
3130 {
3131   // recreate folder avery time and export objects to new one
3132   //
3133   AlidNdPtAnalysispPb * comp=this;
3134   TFolder *folder = comp->GetAnalysisFolder();
3135
3136   TString name, title;
3137   TFolder *newFolder = 0;
3138   Int_t i = 0;
3139   Int_t size = array->GetSize();
3140
3141   if(folder) {
3142      // get name and title from old folder
3143      name = folder->GetName();
3144      title = folder->GetTitle();
3145
3146          // delete old one
3147      delete folder;
3148
3149          // create new one
3150      newFolder = CreateFolder(name.Data(),title.Data());
3151      newFolder->SetOwner();
3152
3153          // add objects to folder
3154      while(i < size) {
3155            newFolder->Add(array->At(i));
3156            i++;
3157          }
3158   }
3159
3160 return newFolder;
3161 }
3162
3163 //_____________________________________________________________________________
3164 TFolder* AlidNdPtAnalysispPb::CreateFolder(TString name,TString title) {
3165 // create folder for analysed histograms
3166 //
3167 TFolder *folder = 0;
3168   folder = new TFolder(name.Data(),title.Data());
3169
3170   return folder;
3171 }