PWG0/dNdPt/AlidNdPtAnalysisPbPb.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 // AlidNdPtAnalysisPbPb 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 // Author: J.Otwinowski 04/11/2008
  28 // last change: 2011-04-04 by M.Knichel
  29 //------------------------------------------------------------------------------
  30
  31 #include "TH1.h"
  32 #include "TH2.h"
  33 #include "TCanvas.h"
  34 #include "THnSparse.h"
  35
  36 #include "AliHeader.h"
  37 #include "AliGenEventHeader.h"
  38 #include "AliInputEventHandler.h"
  39 #include "AliAnalysisManager.h"
  40 #include "AliStack.h"
  41 #include "AliESDEvent.h"
  42 #include "AliMCEvent.h"
  43 #include "AliESDtrackCuts.h"
  44 #include "AliLog.h"
  45 #include "AliMultiplicity.h"
  46 #include "AliTracker.h"
  47
  48 #include "AliCentrality.h"
  49
  50 #include "AlidNdPtEventCuts.h"
  51 #include "AlidNdPtAcceptanceCuts.h"
  52 #include "AliPhysicsSelection.h"
  53 #include "AliTriggerAnalysis.h"
  54
  55 #include "AliPWG0Helper.h"
  56 #include "AlidNdPtHelper.h"
  57 #include "AlidNdPtAnalysisPbPb.h"
  58
  59
  60 using namespace std;
  61
  62 ClassImp(AlidNdPtAnalysisPbPb)
  63
  64 //_____________________________________________________________________________
  65   AlidNdPtAnalysisPbPb::AlidNdPtAnalysisPbPb(): AlidNdPt(),
  66   fAnalysisFolder(0),
  67   fHistogramsOn(kFALSE),
  68
  69   // rec. track pt vs true track pt correlation matrix
  70   fTrackPtCorrelationMatrix(0),
  71
  72   // event level correction
  73   fGenEventMatrix(0),
  74   fTriggerEventMatrix(0),
  75   fRecEventMatrix(0),
  76
  77   //
  78   // track-event level correction
  79   //
  80   fGenTrackEventMatrix(0),
  81
  82   fTriggerTrackEventMatrix(0),
  83
  84   fRecTrackEventMatrix(0),
  85
  86   // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
  87   fGenTrackMatrix(0),
  88   fGenPrimTrackMatrix(0),
  89   fRecPrimTrackMatrix(0),
  90
  91   // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)
  92   fRecTrackMatrix(0),
  93   fRecSecTrackMatrix(0),
  94
  95   // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix)
  96   fRecMultTrackMatrix(0),
  97
  98   // event control histograms
  99   fMCEventHist1(0),
 100   fRecEventHist1(0),
 101   fRecEventHist2(0),
 102   fRecMCEventHist1(0),
 103   fRecMCEventHist2(0),
 104
 105   // rec. pt and eta resolution w.r.t MC
 106   fRecMCTrackHist1(0),
 107
 108   //multple reconstructed tracks
 109   fMCMultRecTrackHist1(0),
 110
 111   // rec. track control histograms
 112   fRecTrackHist3(0),
 113
 114   fTriggerAnalysis(0),
 115
 116   fCentralityEstimator(0),
 117
 118   fMultNbins(0),
 119   fPtNbins(0),
 120   fPtCorrNbins(0),
 121   fEtaNbins(0),
 122   fZvNbins(0),
 123   fCentralityNbins(0),
 124   fBinsMult(0),
 125   fBinsPt(0),
 126   fBinsPtCorr(0),
 127   fBinsEta(0),
 128   fBinsZv(0),
 129   fBinsCentrality(0),
 130
 131   fIsInit(kFALSE)
 132 {
 133   // default constructor
 134   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
 135     fMCTrackHist1[i]=0;
 136     fMCPrimTrackHist1[i]=0;
 137     fMCPrimTrackHist2[i]=0;
 138     fMCSecTrackHist1[i]=0;
 139     fRecTrackHist1[i]=0;
 140     fRecTrackHist2[i]=0;
 141     fRecTrackMultHist1[i]=0;
 142   }
 143   //Init();
 144   SetCentralityEstimator();
 145 }
 146
 147 //_____________________________________________________________________________
 148 AlidNdPtAnalysisPbPb::AlidNdPtAnalysisPbPb(Char_t* name, Char_t* title): AlidNdPt(name,title),
 149   fAnalysisFolder(0),
 150   fHistogramsOn(kFALSE),
 151
 152   // rec. track pt vs true track pt correlation matrix
 153   fTrackPtCorrelationMatrix(0),
 154
 155   // event level correction
 156   fGenEventMatrix(0),
 157   fTriggerEventMatrix(0),
 158   fRecEventMatrix(0),
 159
 160   //
 161   // track-event level correction
 162   //
 163   fGenTrackEventMatrix(0),
 164
 165   fTriggerTrackEventMatrix(0),
 166
 167   fRecTrackEventMatrix(0),
 168
 169   // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
 170   fGenTrackMatrix(0),
 171   fGenPrimTrackMatrix(0),
 172   fRecPrimTrackMatrix(0),
 173
 174   // secondary track contamination correction (fRecTrackMatrix - fRecSecTrackMatrix)
 175   fRecTrackMatrix(0),
 176   fRecSecTrackMatrix(0),
 177
 178   // multiple rec. track contamination corrections (fRecTrackMatrix - fRecMultTrackMatrix)
 179   fRecMultTrackMatrix(0),
 180
 181
 182
 183
 184   // event control histograms
 185   fMCEventHist1(0),
 186   fRecEventHist1(0),
 187   fRecEventHist2(0),
 188   fRecMCEventHist1(0),
 189   fRecMCEventHist2(0),
 190
 191   // rec. pt and eta resolution w.r.t MC
 192   fRecMCTrackHist1(0),
 193
 194   //multiple reconstructed tracks
 195   fMCMultRecTrackHist1(0),
 196
 197   // rec. track control histograms
 198   fRecTrackHist3(0),
 199
 200   fTriggerAnalysis(0),
 201
 202   fCentralityEstimator(0),
 203
 204   fMultNbins(0),
 205   fPtNbins(0),
 206   fPtCorrNbins(0),
 207   fEtaNbins(0),
 208   fZvNbins(0),
 209   fCentralityNbins(0),
 210   fBinsMult(0),
 211   fBinsPt(0),
 212   fBinsPtCorr(0),
 213   fBinsEta(0),
 214   fBinsZv(0),
 215   fBinsCentrality(0),
 216
 217   fIsInit(kFALSE)
 218
 219 {
 220   //
 221   // constructor
 222   //
 223   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
 224     fMCTrackHist1[i]=0;
 225     fMCPrimTrackHist1[i]=0;
 226     fMCPrimTrackHist2[i]=0;
 227     fMCSecTrackHist1[i]=0;
 228     fRecTrackHist1[i]=0;
 229     fRecTrackHist2[i]=0;
 230     fRecTrackMultHist1[i]=0;
 231   }
 232
 233  // Init();
 234  SetCentralityEstimator();
 235 }
 236
 237 //_____________________________________________________________________________
 238 AlidNdPtAnalysisPbPb::~AlidNdPtAnalysisPbPb() {
 239   //
 240   // destructor
 241   //
 242   if(fTrackPtCorrelationMatrix) delete fTrackPtCorrelationMatrix; fTrackPtCorrelationMatrix=0;
 243   //
 244   if(fGenEventMatrix) delete fGenEventMatrix; fGenEventMatrix=0;
 245   if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;
 246   if(fRecEventMatrix) delete fRecEventMatrix; fRecEventMatrix=0;
 247
 248   //
 249   if(fGenTrackEventMatrix) delete fGenTrackEventMatrix; fGenTrackEventMatrix=0;
 250   if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;
 251   if(fRecTrackEventMatrix) delete fRecTrackEventMatrix; fRecTrackEventMatrix=0;
 252
 253   //
 254   if(fGenTrackMatrix) delete fGenTrackMatrix; fGenTrackMatrix=0;
 255   if(fGenPrimTrackMatrix) delete fGenPrimTrackMatrix; fGenPrimTrackMatrix=0;
 256   if(fRecPrimTrackMatrix) delete fRecPrimTrackMatrix; fRecPrimTrackMatrix=0;
 257   //
 258   if(fRecTrackMatrix) delete fRecTrackMatrix; fRecTrackMatrix=0;
 259   if(fRecSecTrackMatrix) delete fRecSecTrackMatrix; fRecSecTrackMatrix=0;
 260   //
 261   if(fRecMultTrackMatrix) delete fRecMultTrackMatrix; fRecMultTrackMatrix=0;
 262   //
 263   // Control histograms
 264   //
 265   if(fMCEventHist1) delete fMCEventHist1; fMCEventHist1=0;
 266   if(fRecEventHist1) delete fRecEventHist1; fRecEventHist1=0;
 267   if(fRecEventHist2) delete fRecEventHist2; fRecEventHist2=0;
 268   if(fRecMCEventHist1) delete fRecMCEventHist1; fRecMCEventHist1=0;
 269   if(fRecMCEventHist2) delete fRecMCEventHist2; fRecMCEventHist2=0;
 270
 271   //
 272   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
 273     if(fMCTrackHist1[i]) delete fMCTrackHist1[i]; fMCTrackHist1[i]=0;
 274     if(fMCPrimTrackHist1[i]) delete fMCPrimTrackHist1[i]; fMCPrimTrackHist1[i]=0;
 275     if(fMCPrimTrackHist2[i]) delete fMCPrimTrackHist2[i]; fMCPrimTrackHist2[i]=0;
 276     if(fMCSecTrackHist1[i]) delete fMCSecTrackHist1[i]; fMCSecTrackHist1[i]=0;
 277     if(fRecTrackHist1[i]) delete fRecTrackHist1[i]; fRecTrackHist1[i]=0;
 278     if(fRecTrackHist2[i]) delete fRecTrackHist2[i]; fRecTrackHist2[i]=0;
 279     if(fRecTrackMultHist1[i]) delete fRecTrackMultHist1[i]; fRecTrackMultHist1[i]=0;
 280   }
 281   if(fRecMCTrackHist1) delete fRecMCTrackHist1; fRecMCTrackHist1=0;
 282   if(fMCMultRecTrackHist1) delete fMCMultRecTrackHist1; fMCMultRecTrackHist1=0;
 283   if(fRecTrackHist3) delete fRecTrackHist3; fRecTrackHist3=0;
 284   //
 285   if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;
 286
 287   if (fTriggerAnalysis) delete fTriggerAnalysis;  fTriggerAnalysis = 0;
 288
 289   if (fBinsMult) delete fBinsMult; fBinsMult=0;
 290   if (fBinsPt) delete fBinsPt; fBinsPt=0;
 291   if (fBinsPtCorr) delete fBinsPtCorr; fBinsPtCorr=0;
 292   if (fBinsEta) delete fBinsEta; fBinsEta=0;
 293   if (fBinsZv) delete fBinsZv; fBinsZv=0;
 294   if (fBinsCentrality) delete fBinsCentrality; fBinsCentrality=0;
 295 }
 296
 297 //_____________________________________________________________________________
 298 void AlidNdPtAnalysisPbPb::Init() {
 299
 300   Int_t multNbins = 47;
 301   Int_t ptNbinsTrackEventCorr = 37;
 302   Int_t ptNbins = 68;
 303   Int_t etaNbins = 30;
 304   Int_t zvNbins = 12;
 305   Int_t centralityBins = 11;
 306
 307   Double_t binsMultDefault[48] = {
 308                    -0.5, 0.5 , 1.5 , 2.5 , 3.5 , 4.5 , 5.5 , 6.5 , 7.5 , 8.5,
 309                     9.5, 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5,
 310                     19.5, 20.5, 30.5, 40.5 , 50.5 , 60.5 , 70.5 , 80.5 , 90.5 , 100.5,
 311                     200.5, 300.5, 400.5, 500.5, 600.5, 700.5, 800.5, 900.5, 1000.5, 2000.5,
 312                     3000.5, 4000.5, 5000.5, 6000.5, 7000.5, 8000.5, 9000.5, 10000.5 }; // forPbPb
 313   Double_t binsPtTrackEventCorrDefault[38] = {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,20.0,50.0};
 314   Double_t binsPtDefault[69] = {
 315         0.0, 0.05, 0.1, 0.15, 0.2, 0.25, 0.3, 0.35, 0.4, 0.45,
 316         0.5, 0.55, 0.6, 0.65, 0.7, 0.75, 0.8, 0.85, 0.9, 0.95,
 317         1.0, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9,
 318         2.0, 2.2, 2.4, 2.6, 2.8, 3.0, 3.2, 3.4, 3.6, 3.8,
 319         4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.0, 8.0, 9.0, 10.0,
 320         11.0, 12.0, 13.0, 14.0, 15.0, 16.0, 18.0, 20.0, 22.0, 24.0,
 321         26.0, 28.0, 30.0, 32.0, 34.0, 36.0, 40.0, 45.0, 50.0 };
 322
 323   Double_t binsEtaDefault[31] = {
 324         -1.5, -1.4, -1.3, -1.2, -1.1, -1.0, -0.9, -0.8, -0.7, -0.6,
 325         -0.5, -0.4, -0.3, -0.2, -0.1,  0.0,  0.1,  0.2,  0.3,  0.4,
 326          0.5,  0.6,  0.7,  0.8,  0.9,  1.0,  1.1,  1.2,  1.3,  1.4,
 327          1.5};
 328
 329   Double_t binsZvDefault[13] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
 330
 331   Double_t binsCentralityDefault[12] = {0., 5., 10., 20., 30., 40., 50., 60., 70., 80., 90., 100.};
 332
 333   Double_t* binsMult = binsMultDefault;
 334   Double_t* binsPtTrackEventCorr = binsPtTrackEventCorrDefault;
 335   Double_t* binsPt = binsPtDefault;
 336   Double_t* binsEta = binsEtaDefault;
 337   Double_t* binsZv = binsZvDefault;
 338   Double_t* binsCentrality = binsCentralityDefault;
 339
 340   if (fMultNbins > 0) { multNbins = fMultNbins; binsMult = fBinsMult; }
 341   if (fPtNbins > 0) { ptNbins  = fPtNbins; binsPt = fBinsPt; }
 342   if (fPtCorrNbins > 0) { ptNbinsTrackEventCorr  = fPtCorrNbins; binsPtTrackEventCorr = fBinsPtCorr;  }
 343   if (fEtaNbins > 0) { etaNbins  = fEtaNbins; binsEta = fBinsEta; }
 344   if (fZvNbins > 0) { zvNbins  = fZvNbins; binsZv = fBinsZv; }
 345   if (fCentralityNbins > 0) { centralityBins  = fCentralityNbins; binsCentrality = fBinsCentrality; }
 346
 347   Int_t binsTrackEventCorrMatrix[4]={zvNbins,ptNbinsTrackEventCorr,etaNbins,centralityBins};
 348   Int_t binsTrackEvent[4]={zvNbins,ptNbins,etaNbins,centralityBins};
 349   Int_t binsTrackPtCorrelationMatrix[4]={ptNbins,ptNbins,etaNbins,centralityBins};
 350
 351
 352   fTrackPtCorrelationMatrix = new THnSparseF("fTrackPtCorrelationMatrix","Pt:mcPt:mcEta:Centrality",4,binsTrackPtCorrelationMatrix);
 353   fTrackPtCorrelationMatrix->SetBinEdges(0,binsPt);
 354   fTrackPtCorrelationMatrix->SetBinEdges(1,binsPt);
 355   fTrackPtCorrelationMatrix->SetBinEdges(2,binsEta);
 356   fTrackPtCorrelationMatrix->SetBinEdges(3,binsCentrality);
 357   fTrackPtCorrelationMatrix->GetAxis(0)->SetTitle("Pt (GeV/c)");
 358   fTrackPtCorrelationMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 359   fTrackPtCorrelationMatrix->GetAxis(2)->SetTitle("mcEta");
 360   fTrackPtCorrelationMatrix->GetAxis(3)->SetTitle("Centrality");
 361   fTrackPtCorrelationMatrix->Sumw2();
 362
 363   //
 364   // Efficiency and contamination correction matrices
 365   //
 366   Int_t binsEventMatrix[3]={zvNbins,multNbins,centralityBins};
 367   Double_t minEventMatrix[3]={-30.,-0.5,0.};
 368   Double_t maxEventMatrix[3]={30.,10000.5,100.};
 369
 370   fGenEventMatrix = new THnSparseF("fGenEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix,minEventMatrix,maxEventMatrix);
 371   fGenEventMatrix->SetBinEdges(0,binsZv);
 372   fGenEventMatrix->SetBinEdges(1,binsMult);
 373   fGenEventMatrix->SetBinEdges(2,binsCentrality);
 374   fGenEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 375   fGenEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");
 376   fGenEventMatrix->GetAxis(2)->SetTitle("Centrality");
 377   fGenEventMatrix->Sumw2();
 378   //
 379   fTriggerEventMatrix = new THnSparseF("fTriggerEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix,minEventMatrix,maxEventMatrix);
 380   fTriggerEventMatrix->SetBinEdges(0,binsZv);
 381   fTriggerEventMatrix->SetBinEdges(1,binsMult);
 382   fTriggerEventMatrix->SetBinEdges(2,binsCentrality);
 383   fTriggerEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 384   fTriggerEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");
 385   fTriggerEventMatrix->GetAxis(2)->SetTitle("Centrality");
 386   fTriggerEventMatrix->Sumw2();
 387   //
 388   fRecEventMatrix = new THnSparseF("fRecEventMatrix","mcZv:multMB:Centrality",3,binsEventMatrix,minEventMatrix,maxEventMatrix);
 389   fRecEventMatrix->SetBinEdges(0,binsZv);
 390   fRecEventMatrix->SetBinEdges(1,binsMult);
 391   fRecEventMatrix->SetBinEdges(2,binsCentrality);
 392   fRecEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 393   fRecEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");
 394   fRecEventMatrix->GetAxis(2)->SetTitle("Centrality");
 395   fRecEventMatrix->Sumw2();
 396
 397   //
 398   // track to event corrections
 399   //
 400   fGenTrackEventMatrix = new THnSparseF("fGenTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 401   fGenTrackEventMatrix->SetBinEdges(0,binsZv);
 402   fGenTrackEventMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 403   fGenTrackEventMatrix->SetBinEdges(2,binsEta);
 404   fGenTrackEventMatrix->SetBinEdges(3,binsCentrality);
 405   fGenTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 406   fGenTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 407   fGenTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
 408   fGenTrackEventMatrix->GetAxis(3)->SetTitle("Centrality");
 409   fGenTrackEventMatrix->Sumw2();
 410   //
 411   fTriggerTrackEventMatrix = new THnSparseF("fTriggerTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 412   fTriggerTrackEventMatrix->SetBinEdges(0,binsZv);
 413   fTriggerTrackEventMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 414   fTriggerTrackEventMatrix->SetBinEdges(2,binsEta);
 415   fTriggerTrackEventMatrix->SetBinEdges(3,binsCentrality);
 416   fTriggerTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 417   fTriggerTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 418   fTriggerTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
 419   fTriggerTrackEventMatrix->GetAxis(3)->SetTitle("Centrality");
 420   fTriggerTrackEventMatrix->Sumw2();
 421   //
 422   fRecTrackEventMatrix = new THnSparseF("fRecTrackEventMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 423   fRecTrackEventMatrix->SetBinEdges(0,binsZv);
 424   fRecTrackEventMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 425   fRecTrackEventMatrix->SetBinEdges(2,binsEta);
 426   fRecTrackEventMatrix->SetBinEdges(3,binsCentrality);
 427   fRecTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 428   fRecTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 429   fRecTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");
 430   fRecTrackEventMatrix->GetAxis(3)->SetTitle("Centrality");
 431   fRecTrackEventMatrix->Sumw2();
 432
 433   //
 434   // tracks correction matrices
 435   //
 436   fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 437   fGenTrackMatrix->SetBinEdges(0,binsZv);
 438   fGenTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 439   fGenTrackMatrix->SetBinEdges(2,binsEta);
 440   fGenTrackMatrix->SetBinEdges(3,binsCentrality);
 441   fGenTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 442   fGenTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 443   fGenTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 444   fGenTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 445   fGenTrackMatrix->Sumw2();
 446
 447   fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 448   fGenPrimTrackMatrix->SetBinEdges(0,binsZv);
 449   fGenPrimTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 450   fGenPrimTrackMatrix->SetBinEdges(2,binsEta);
 451   fGenPrimTrackMatrix->SetBinEdges(3,binsCentrality);
 452   fGenPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 453   fGenPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 454   fGenPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 455   fGenPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 456   fGenPrimTrackMatrix->Sumw2();
 457
 458
 459   fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 460   fRecPrimTrackMatrix->SetBinEdges(0,binsZv);
 461   fRecPrimTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 462   fRecPrimTrackMatrix->SetBinEdges(2,binsEta);
 463   fRecPrimTrackMatrix->SetBinEdges(3,binsCentrality);
 464   fRecPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 465   fRecPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 466   fRecPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 467   fRecPrimTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 468   fRecPrimTrackMatrix->Sumw2();
 469
 470   //
 471   fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 472   fRecTrackMatrix->SetBinEdges(0,binsZv);
 473   fRecTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 474   fRecTrackMatrix->SetBinEdges(2,binsEta);
 475   fRecTrackMatrix->SetBinEdges(3,binsCentrality);
 476   fRecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 477   fRecTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 478   fRecTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 479   fRecTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 480   fRecTrackMatrix->Sumw2();
 481
 482   fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 483   fRecSecTrackMatrix->SetBinEdges(0,binsZv);
 484   fRecSecTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 485   fRecSecTrackMatrix->SetBinEdges(2,binsEta);
 486   fRecSecTrackMatrix->SetBinEdges(3,binsCentrality);
 487   fRecSecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 488   fRecSecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)");
 489   fRecSecTrackMatrix->GetAxis(2)->SetTitle("Eta");
 490   fRecSecTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 491   fRecSecTrackMatrix->Sumw2();
 492
 493   //
 494   fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta:Centrality",4,binsTrackEventCorrMatrix);
 495   fRecMultTrackMatrix->SetBinEdges(0,binsZv);
 496   fRecMultTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);
 497   fRecMultTrackMatrix->SetBinEdges(2,binsEta);
 498   fRecMultTrackMatrix->SetBinEdges(3,binsCentrality);
 499   fRecMultTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");
 500   fRecMultTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");
 501   fRecMultTrackMatrix->GetAxis(2)->SetTitle("mcEta");
 502   fRecMultTrackMatrix->GetAxis(3)->SetTitle("Centrality");
 503   fRecMultTrackMatrix->Sumw2();
 504
 505   //
 506   // Control analysis histograms
 507   //
 508   Int_t binsMCEventHist1[4]={100,100,zvNbins,centralityBins};
 509   Double_t minMCEventHist1[4]={-0.1,-0.1,-30.,0.};
 510   Double_t maxMCEventHist1[4]={0.1,0.1,30.,100.};
 511   fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv:Centrality",4,binsMCEventHist1,minMCEventHist1,maxMCEventHist1);
 512   fMCEventHist1->SetBinEdges(2,binsZv);
 513   fMCEventHist1->SetBinEdges(3,binsCentrality);
 514   fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)");
 515   fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)");
 516   fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)");
 517   fMCEventHist1->GetAxis(3)->SetTitle("Centrality");
 518   fMCEventHist1->Sumw2();
 519
 520   //
 521   Int_t binsRecEventHist1[4]={100,100,zvNbins,centralityBins};
 522   Double_t minRecEventHist1[4]={-3.,-3.,-30.,0.};
 523   Double_t maxRecEventHist1[4]={3.,3.,30.,100.};
 524
 525   fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv:Centrality",4,binsRecEventHist1,minRecEventHist1,maxRecEventHist1);
 526   fRecEventHist1->SetBinEdges(2,binsZv);
 527   fRecEventHist1->SetBinEdges(3,binsCentrality);
 528   fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)");
 529   fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)");
 530   fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)");
 531   fRecEventHist1->GetAxis(3)->SetTitle("Centrality");
 532   fRecEventHist1->Sumw2();
 533
 534   //
 535   Int_t binsRecEventHist2[3]={zvNbins, multNbins, centralityBins};
 536   Double_t minRecEventHist2[3]={-30., -0.5, 0.};
 537   Double_t maxRecEventHist2[3]={30., 10000.5, 100.};
 538
 539   fRecEventHist2 = new THnSparseF("fRecEventHist2","Zv:multMB:Centrality",3,binsRecEventHist2,minRecEventHist2,maxRecEventHist2);
 540   fRecEventHist2->SetBinEdges(0,binsZv);
 541   fRecEventHist2->SetBinEdges(1,binsMult);
 542   fRecEventHist2->SetBinEdges(2,binsCentrality);
 543   fRecEventHist2->GetAxis(0)->SetTitle("Zv (cm)");
 544   fRecEventHist2->GetAxis(1)->SetTitle("multiplicity MB");
 545   fRecEventHist2->GetAxis(2)->SetTitle("Centrality");
 546   fRecEventHist2->Sumw2();
 547
 548   //
 549   Double_t kFact = 0.1;
 550   Int_t binsRecMCEventHist1[4]={100,100,100, centralityBins};
 551   Double_t minRecMCEventHist1[4]={-10.0*kFact,-10.0*kFact,-10.0*kFact, 0.};
 552   Double_t maxRecMCEventHist1[4]={10.0*kFact,10.0*kFact,10.0*kFact, 100.};
 553
 554   fRecMCEventHist1 = new THnSparseF("fRecMCEventHist1","Xv-mcXv:Yv-mcYv:Zv-mcZv:Centrality",4,binsRecMCEventHist1,minRecMCEventHist1,maxRecMCEventHist1);
 555   fRecMCEventHist1->SetBinEdges(3,binsCentrality);
 556   fRecMCEventHist1->GetAxis(0)->SetTitle("Xv-mcXv (cm)");
 557   fRecMCEventHist1->GetAxis(1)->SetTitle("Yv-mcYv (cm)");
 558   fRecMCEventHist1->GetAxis(2)->SetTitle("Zv-mcZv (cm)");
 559   fRecMCEventHist1->GetAxis(3)->SetTitle("Centrality");
 560
 561   fRecMCEventHist1->Sumw2();
 562
 563   //
 564   Int_t binsRecMCEventHist2[4]={100,100,multNbins, centralityBins};
 565   Double_t minRecMCEventHist2[4]={-10.0*kFact,-10.0*kFact,-0.5, 0.};
 566   Double_t maxRecMCEventHist2[4]={10.0*kFact,10.0*kFact,10000.5, 100.};
 567
 568   fRecMCEventHist2 = new THnSparseF("fRecMCEventHist2","Xv-mcXv:Zv-mcZv:mult:Centrality",4,binsRecMCEventHist2,minRecMCEventHist2,maxRecMCEventHist2);
 569   fRecMCEventHist2->SetBinEdges(2,binsMult);
 570   fRecMCEventHist2->SetBinEdges(3,binsCentrality);
 571   fRecMCEventHist2->GetAxis(0)->SetTitle("Xv-mcXv (cm)");
 572   fRecMCEventHist2->GetAxis(1)->SetTitle("Zv-mcZv (cm)");
 573   fRecMCEventHist2->GetAxis(2)->SetTitle("multiplicity");
 574   fRecMCEventHist2->GetAxis(3)->SetTitle("Centrality");
 575   fRecMCEventHist2->Sumw2();
 576   //
 577   char name[256];
 578   char title[256];
 579   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++)
 580   {
 581   // THnSparse track histograms
 582
 583   Int_t binsMCTrackHist1[4]=  {ptNbins, etaNbins, 90, centralityBins};
 584   Double_t minMCTrackHist1[4]={0.,-1.5,0., 0.};
 585   Double_t maxMCTrackHist1[4]={50,1.5,2.*TMath::Pi(), 100.};
 586   sprintf(name,"fMCTrackHist1_%d",i);
 587   sprintf(title,"mcPt:mcEta:mcPhi:Centrality");
 588
 589   fMCTrackHist1[i] = new THnSparseF(name,title,4,binsMCTrackHist1,minMCTrackHist1,maxMCTrackHist1);
 590   fMCTrackHist1[i]->SetBinEdges(0,binsPt);
 591   fMCTrackHist1[i]->SetBinEdges(1,binsEta);
 592   fMCTrackHist1[i]->SetBinEdges(3,binsCentrality);
 593   fMCTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
 594   fMCTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
 595   fMCTrackHist1[i]->GetAxis(2)->SetTitle("mcPhi (rad)");
 596   fMCTrackHist1[i]->GetAxis(3)->SetTitle("Centrality");
 597   fMCTrackHist1[i]->Sumw2();
 598
 599   Int_t binsMCPrimTrackHist1[6]=  {ptNbins,etaNbins,6,20,4000, centralityBins};
 600   Double_t minMCPrimTrackHist1[6]={0.,-1.5,0.,0.,0., 0.};
 601   Double_t maxMCPrimTrackHist1[6]={50.,1.5,6.,20.,4000., 100.};
 602   sprintf(name,"fMCPrimTrackHist1_%d",i);
 603   sprintf(title,"mcPt:mcEta:pid:mech:mother:Centrality");
 604
 605   fMCPrimTrackHist1[i] = new THnSparseF(name,title,6,binsMCPrimTrackHist1,minMCPrimTrackHist1,maxMCPrimTrackHist1);
 606   fMCPrimTrackHist1[i]->SetBinEdges(0,binsPt);
 607   fMCPrimTrackHist1[i]->SetBinEdges(1,binsEta);
 608   fMCPrimTrackHist1[i]->SetBinEdges(5,binsCentrality);
 609   fMCPrimTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
 610   fMCPrimTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
 611   fMCPrimTrackHist1[i]->GetAxis(2)->SetTitle("pid");
 612   fMCPrimTrackHist1[i]->GetAxis(3)->SetTitle("mech");
 613   fMCPrimTrackHist1[i]->GetAxis(4)->SetTitle("mother");
 614   fMCPrimTrackHist1[i]->GetAxis(5)->SetTitle("Centrality");
 615   fMCPrimTrackHist1[i]->Sumw2();
 616
 617   Int_t binsMCPrimTrackHist2[4]=  {4000,20,4000,centralityBins};
 618   Double_t minMCPrimTrackHist2[4]={0.,0.,0., 0.};
 619   Double_t maxMCPrimTrackHist2[4]={4000.,20.,4000., 100.};
 620   sprintf(name,"fMCPrimTrackHist2_%d",i);
 621   sprintf(title,"pdg:mech:mother:Centrality");
 622
 623   fMCPrimTrackHist2[i] = new THnSparseF(name,title,4,binsMCPrimTrackHist2,minMCPrimTrackHist2,maxMCPrimTrackHist2);
 624   fMCPrimTrackHist2[i]->SetBinEdges(3,binsCentrality);
 625   fMCPrimTrackHist2[i]->GetAxis(0)->SetTitle("pdg");
 626   fMCPrimTrackHist2[i]->GetAxis(1)->SetTitle("mech");
 627   fMCPrimTrackHist2[i]->GetAxis(2)->SetTitle("mother");
 628   fMCPrimTrackHist2[i]->GetAxis(3)->SetTitle("Centrality");
 629   fMCPrimTrackHist2[i]->Sumw2();
 630
 631   Int_t binsMCSecTrackHist1[6]=  {ptNbins,etaNbins,6,20,4000, centralityBins};
 632   Double_t minMCSecTrackHist1[6]={0.,-1.5,0.,0.,0., 0.};
 633   Double_t maxMCSecTrackHist1[6]={50.,1.5,6.,20.,4000., 100.};
 634   sprintf(name,"fMCSecTrackHist1_%d",i);
 635   sprintf(title,"mcPt:mcEta:pid:mech:mother:Centrality");
 636
 637   fMCSecTrackHist1[i] = new THnSparseF(name,title,6,binsMCSecTrackHist1,minMCSecTrackHist1,maxMCSecTrackHist1);
 638   fMCSecTrackHist1[i]->SetBinEdges(0,binsPt);
 639   fMCSecTrackHist1[i]->SetBinEdges(1,binsEta);
 640   fMCSecTrackHist1[i]->SetBinEdges(5,binsCentrality);
 641   fMCSecTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");
 642   fMCSecTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");
 643   fMCSecTrackHist1[i]->GetAxis(2)->SetTitle("pid");
 644   fMCSecTrackHist1[i]->GetAxis(3)->SetTitle("mech");
 645   fMCSecTrackHist1[i]->GetAxis(4)->SetTitle("mother");
 646   fMCSecTrackHist1[i]->GetAxis(5)->SetTitle("Centrality");
 647   fMCSecTrackHist1[i]->Sumw2();
 648
 649   Int_t binsRecTrackHist1[4]={ptNbins,etaNbins,90, centralityBins};
 650   Double_t minRecTrackHist1[4]={0.,-1.5,0., 0.};
 651   Double_t maxRecTrackHist1[4]={50.,1.5,2.*TMath::Pi(), 100.};
 652   sprintf(name,"fRecTrackHist1_%d",i);
 653   sprintf(title,"Pt:Eta:Phi:Centrality");
 654   fRecTrackHist1[i] = new THnSparseF(name,title,4,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1);
 655   fRecTrackHist1[i]->SetBinEdges(0,binsPt);
 656   fRecTrackHist1[i]->SetBinEdges(1,binsEta);
 657   fRecTrackHist1[i]->SetBinEdges(3,binsCentrality);
 658   fRecTrackHist1[i]->GetAxis(0)->SetTitle("p_{T} (GeV/c)");
 659   fRecTrackHist1[i]->GetAxis(1)->SetTitle("#eta");
 660   fRecTrackHist1[i]->GetAxis(2)->SetTitle("#phi (rad)");
 661   fRecTrackHist1[i]->GetAxis(3)->SetTitle("Centrality");
 662   fRecTrackHist1[i]->Sumw2();
 663
 664   sprintf(name,"fRecTrackHist2_%d",i);
 665   sprintf(title,"Zv:Pt:Eta:Centrality");
 666   fRecTrackHist2[i] = new THnSparseF(name,title,4,binsTrackEvent);
 667   fRecTrackHist2[i]->SetBinEdges(0,binsZv);
 668   fRecTrackHist2[i]->SetBinEdges(1,binsPt);
 669   fRecTrackHist2[i]->SetBinEdges(2,binsEta);
 670   fRecTrackHist2[i]->SetBinEdges(3,binsCentrality);
 671   fRecTrackHist2[i]->GetAxis(0)->SetTitle("Zv (cm)");
 672   fRecTrackHist2[i]->GetAxis(1)->SetTitle("p_{T} (GeV/c)");
 673   fRecTrackHist2[i]->GetAxis(2)->SetTitle("#eta");
 674   fRecTrackHist2[i]->GetAxis(3)->SetTitle("Centrality");
 675   fRecTrackHist2[i]->Sumw2();
 676
 677   //
 678   Int_t binsRecTrackMultHist1[3]={ptNbins,multNbins, centralityBins};
 679   Double_t minRecTrackMultHist1[3]={0.,-0.5, -0.};
 680   Double_t maxRecTrackMultHist1[3]={50.,10000.5, 100.};
 681   sprintf(name,"fRecTrackMultHist_%d",i);
 682   sprintf(title,"Pt:Mult:Centrality");
 683   fRecTrackMultHist1[i] = new THnSparseF(name,title,3,binsRecTrackMultHist1,minRecTrackMultHist1,maxRecTrackMultHist1);
 684   fRecTrackMultHist1[i]->SetBinEdges(0,binsPt);
 685   fRecTrackMultHist1[i]->SetBinEdges(1,binsMult);
 686   fRecTrackMultHist1[i]->SetBinEdges(2,binsCentrality);
 687   fRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
 688   fRecTrackMultHist1[i]->GetAxis(1)->SetTitle("multiplicity");
 689   fRecTrackMultHist1[i]->GetAxis(2)->SetTitle("Centrality");
 690   fRecTrackMultHist1[i]->Sumw2();
 691   }
 692
 693   Int_t binsRecMCTrackHist1[5] = {ptNbins,etaNbins,100,100, centralityBins};
 694   Double_t minRecMCTrackHist1[5]={0.,-1.5,-0.5,-0.5, 0.};
 695   Double_t maxRecMCTrackHist1[5]={50.,1.5,0.5,0.5, 100.};
 696
 697   sprintf(name,"fRecMCTrackHist1");
 698   sprintf(title,"mcPt:mcEta:(Pt-mcPt)/mcPt:(Eta-mcEta):Centrality");
 699   fRecMCTrackHist1 = new THnSparseF(name,title,5,binsRecMCTrackHist1,minRecMCTrackHist1,maxRecMCTrackHist1);
 700   fRecMCTrackHist1->SetBinEdges(0,binsPt);
 701   fRecMCTrackHist1->SetBinEdges(1,binsEta);
 702   fRecMCTrackHist1->SetBinEdges(4,binsCentrality);
 703   fRecMCTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
 704   fRecMCTrackHist1->GetAxis(1)->SetTitle("mcEta");
 705   fRecMCTrackHist1->GetAxis(2)->SetTitle("(Pt-mcPt)/mcPt");
 706   fRecMCTrackHist1->GetAxis(3)->SetTitle("Eta-mcEta");
 707   fRecMCTrackHist1->GetAxis(4)->SetTitle("Centrality");
 708
 709   Int_t binsMCMultRecTrackHist1[4] = {ptNbins,etaNbins,6, centralityBins};
 710   Double_t minMCMultRecTrackHist1[4]={0.,-1.5,0., 0.};
 711   Double_t maxMCMultRecTrackHist1[4]={50.,1.5,6., 100.};
 712   sprintf(name,"fMCMultRecTrackHist1");
 713   sprintf(title,"mcPt:mcEta:pid:Centrality");
 714   fMCMultRecTrackHist1 = new THnSparseF(name,title,4,binsMCMultRecTrackHist1,minMCMultRecTrackHist1,maxMCMultRecTrackHist1);
 715   fMCMultRecTrackHist1->SetBinEdges(0,binsPt);
 716   fMCMultRecTrackHist1->SetBinEdges(1,binsEta);
 717   fMCMultRecTrackHist1->SetBinEdges(3,binsCentrality);
 718   fMCMultRecTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
 719   fMCMultRecTrackHist1->GetAxis(1)->SetTitle("mcEta");
 720   fMCMultRecTrackHist1->GetAxis(2)->SetTitle("pid");
 721   fMCMultRecTrackHist1->GetAxis(3)->SetTitle("Centrality");
 722
 723   //nClust:chi2PerClust:pt:eta:phi:Centrality
 724   Int_t binsRecTrackHist3[6]={160,100,ptNbins,etaNbins,90, centralityBins};
 725   Double_t minRecTrackHist3[6]={0., 0., 0., -1.5, 0., 0.};
 726   Double_t maxRecRecTrackHist3[6]={160.,10., 50., 1.5, 2.*TMath::Pi(), 100.};
 727
 728   fRecTrackHist3 = new THnSparseF("fRecTrackHist3","nClust:chi2PerClust:pt:eta:phi:Centrality",6,binsRecTrackHist3,minRecTrackHist3,maxRecRecTrackHist3);
 729   fRecTrackHist3->SetBinEdges(2,binsPt);
 730   fRecTrackHist3->SetBinEdges(3,binsEta);
 731   fRecTrackHist3->SetBinEdges(5,binsCentrality);
 732   fRecTrackHist3->GetAxis(0)->SetTitle("nClust");
 733   fRecTrackHist3->GetAxis(1)->SetTitle("chi2PerClust");
 734   fRecTrackHist3->GetAxis(2)->SetTitle("p_{T} (GeV/c)");
 735   fRecTrackHist3->GetAxis(3)->SetTitle("#eta");
 736   fRecTrackHist3->GetAxis(4)->SetTitle("#phi (rad)");
 737   fRecTrackHist3->GetAxis(5)->SetTitle("Centrality");
 738   fRecTrackHist3->Sumw2();
 739
 740
 741   // init folder
 742   fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");
 743
 744   // init trigger analysis (for zdc cut)
 745   fTriggerAnalysis = new AliTriggerAnalysis;
 746
 747   // set init flag
 748   fIsInit = kTRUE;
 749 }
 750
 751 //_____________________________________________________________________________
 752 void AlidNdPtAnalysisPbPb::Process(AliESDEvent *const esdEvent, AliMCEvent *const mcEvent)
 753 {
 754   //  init if not done already
 755   if (!fIsInit) { Init(); }
 756   //
 757   // Process real and/or simulated events
 758   //
 759   if(!esdEvent) {
 760     AliDebug(AliLog::kError, "esdEvent not available");
 761     return;
 762   }
 763
 764   // get selection cuts
 765   AlidNdPtEventCuts *evtCuts = GetEventCuts();
 766   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts();
 767   AlidNdPtAcceptanceCuts *recCuts = GetRecAcceptanceCuts();
 768   AliESDtrackCuts *esdTrackCuts = GetTrackCuts();
 769
 770   if(!evtCuts || !accCuts  || !esdTrackCuts) {
 771     AliDebug(AliLog::kError, "cuts not available");
 772     return;
 773   }
 774   if (0 == recCuts) { recCuts = accCuts;}
 775
 776   // trigger selection
 777   Bool_t isEventTriggered = kTRUE;
 778   AliPhysicsSelection *physicsSelection = NULL;
 779   AliTriggerAnalysis* triggerAnalysis = NULL;
 780
 781   //
 782   AliInputEventHandler* inputHandler = (AliInputEventHandler*) AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler();
 783   if (!inputHandler)
 784   {
 785     Printf("ERROR: Could not receive input handler");
 786     return;
 787   }
 788
 789   if(evtCuts->IsTriggerRequired())
 790   {
 791     // always MB
 792     isEventTriggered = inputHandler->IsEventSelected() & AliVEvent::kMB;
 793
 794     physicsSelection = static_cast<AliPhysicsSelection*> (inputHandler->GetEventSelection());
 795     if(!physicsSelection) return;
 796     //SetPhysicsTriggerSelection(physicsSelection);
 797
 798     if (isEventTriggered && (GetTrigger() == AliTriggerAnalysis::kV0AND)) {
 799       // set trigger (V0AND)
 800       triggerAnalysis = physicsSelection->GetTriggerAnalysis();
 801       if(!triggerAnalysis) return;
 802       isEventTriggered = triggerAnalysis->IsOfflineTriggerFired(esdEvent, GetTrigger());
 803     }
 804   }
 805
 806
 807   // centrality determination
 808   Float_t centralityF = -1.;
 809   AliCentrality *esdCentrality = esdEvent->GetCentrality();
 810   centralityF = esdCentrality->GetCentralityPercentile(fCentralityEstimator.Data());
 811
 812   // use MC information
 813   AliHeader* header = 0;
 814   AliGenEventHeader* genHeader = 0;
 815   AliStack* stack = 0;
 816   TArrayF vtxMC(3);
 817
 818   Int_t multMCTrueTracks = 0;
 819   if(IsUseMCInfo())
 820   {
 821     //
 822     if(!mcEvent) {
 823       AliDebug(AliLog::kError, "mcEvent not available");
 824       return;
 825     }
 826     // get MC event header
 827     header = mcEvent->Header();
 828     if (!header) {
 829       AliDebug(AliLog::kError, "Header not available");
 830       return;
 831     }
 832     // MC particle stack
 833     stack = mcEvent->Stack();
 834     if (!stack) {
 835       AliDebug(AliLog::kError, "Stack not available");
 836       return;
 837     }
 838
 839     // get MC vertex
 840     genHeader = header->GenEventHeader();
 841     if (!genHeader) {
 842       AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
 843       return;
 844     }
 845     genHeader->PrimaryVertex(vtxMC);
 846
 847     Double_t vMCEventHist1[4]={vtxMC[0],vtxMC[1],vtxMC[2],centralityF};
 848     fMCEventHist1->Fill(vMCEventHist1);
 849
 850     // multipliticy of all MC primary tracks
 851     // in Zv, pt and eta ranges)
 852     multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);
 853
 854   } // end bUseMC
 855
 856   // get reconstructed vertex
 857   const AliESDVertex* vtxESD = 0;
 858   if(evtCuts->IsRecVertexRequired())
 859   {
 860     if(GetAnalysisMode() == AlidNdPtHelper::kTPC) {
 861       vtxESD = esdEvent->GetPrimaryVertexTPC();
 862     }
 863     else if(GetAnalysisMode() == AlidNdPtHelper::kTPCITS) {
 864       vtxESD = esdEvent->GetPrimaryVertexTracks();
 865     }
 866     else {
 867         return;
 868     }
 869   }
 870
 871   Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD);
 872   //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);
 873   //printf("GetAnalysisMode() %d \n",GetAnalysisMode());
 874
 875   // vertex contributors
 876   Int_t multMBTracks = 0;
 877   if(GetAnalysisMode() == AlidNdPtHelper::kTPC || GetAnalysisMode() == AlidNdPtHelper::kTPCITS)
 878   {
 879      if(vtxESD->GetStatus()) {
 880          multMBTracks = vtxESD->GetNContributors();
 881      }
 882   }
 883   else {
 884     AliDebug(AliLog::kError, Form("Found analysis type %d", GetAnalysisMode()));
 885     return;
 886   }
 887
 888   TObjArray *allChargedTracks=0;
 889   //Int_t multAll=0, multAcc=0, multRec=0;
 890   Int_t multAll=0, multRec=0;
 891   Int_t *labelsAll=0, *labelsAcc=0, *labelsRec=0;
 892
 893   // check event cuts
 894   if(isEventOK && isEventTriggered)
 895   {
 896
 897     // get all charged tracks
 898     allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());
 899     if(!allChargedTracks) return;
 900
 901     Int_t entries = allChargedTracks->GetEntries();
 902
 903     labelsAll = new Int_t[entries];
 904     labelsAcc = new Int_t[entries];
 905     labelsRec = new Int_t[entries];
 906     for(Int_t i=0; i<entries;++i)
 907     {
 908       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
 909
 910       if(!track) continue;
 911       if(track->Charge()==0) continue;
 912
 913
 914       // only postive charged
 915       if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0)
 916         continue;
 917
 918       // only negative charged
 919       if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0)
 920         continue;
 921
 922       //
 923       Double_t values[4] = {vtxESD->GetZv(),track->Pt(),track->Eta(), centralityF};
 924
 925       fRecTrackHist2[AlidNdPtHelper::kAllTracks]->Fill(values);
 926       FillHistograms(track,stack,AlidNdPtHelper::kAllTracks,centralityF);
 927       labelsAll[multAll] = TMath::Abs(track->GetLabel());
 928
 929       multAll++;
 930       if(esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track) && recCuts->AcceptTrackLocalTPC(track)) {
 931
 932          fRecTrackHist2[AlidNdPtHelper::kRecTracks]->Fill(values);
 933          FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,centralityF);
 934          labelsRec[multRec] = TMath::Abs(track->GetLabel());
 935
 936          multRec++;
 937
 938       }
 939      }//loop over entries
 940
 941
 942      // fill track multiplicity histograms
 943      //FillHistograms(allChargedTracks,labelsAll,multAll,labelsAcc,multAcc,labelsRec,multRec);
 944
 945      Double_t vRecEventHist1[4] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv(),centralityF};
 946      fRecEventHist1->Fill(vRecEventHist1);
 947
 948      Double_t vRecEventHist2[3] = {vtxESD->GetZv(),multMBTracks,centralityF};
 949      fRecEventHist2->Fill(vRecEventHist2);
 950
 951    } // triggered and event vertex
 952
 953    if(IsUseMCInfo())
 954    {
 955      //
 956      // event level corrections (zv,N_MB)
 957      //
 958      // all inelastic
 959      Double_t vEventMatrix[3] = {vtxMC[2],multMBTracks,centralityF};
 960      fGenEventMatrix->Fill(vEventMatrix);
 961      if(isEventTriggered) fTriggerEventMatrix->Fill(vEventMatrix);
 962      if(isEventOK && isEventTriggered) fRecEventMatrix->Fill(vEventMatrix);
 963
 964      //
 965      // track-event level corrections (zv,pt,eta)
 966      //
 967      for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
 968      {
 969        TParticle* particle = stack->Particle(iMc);
 970        if (!particle)
 971        continue;
 972
 973        // only charged particles
 974        if(!particle->GetPDG()) continue;
 975        Double_t charge = particle->GetPDG()->Charge()/3.;
 976        if ( TMath::Abs(charge) < 0.001 )
 977         continue;
 978
 979        // only postive charged
 980        if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.)
 981         continue;
 982
 983        // only negative charged
 984        if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.)
 985        continue;
 986
 987        // physical primary
 988        Bool_t prim = stack->IsPhysicalPrimary(iMc);
 989        if(!prim) continue;
 990
 991        // checked accepted
 992        if(accCuts->AcceptTrack(particle))
 993        {
 994         Double_t vTrackEventMatrix[4] = {vtxMC[2], particle->Pt(), particle->Eta(), centralityF};
 995         fGenTrackEventMatrix->Fill(vTrackEventMatrix);
 996
 997         if(!isEventTriggered) continue;
 998         fTriggerTrackEventMatrix->Fill(vTrackEventMatrix);
 999
1000         if(!isEventOK) continue;
1001         fRecTrackEventMatrix->Fill(vTrackEventMatrix);
1002
1003        }// if(accCuts->AcceptTrack(particle))
1004      }// for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
1005
1006      //
1007      // track-level corrections (zv,pt,eta)
1008      //
1009      if(isEventOK && isEventTriggered)
1010      {
1011
1012        // fill MC and rec event control histograms
1013        if(fHistogramsOn) {
1014          Double_t vRecMCEventHist1[4] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2], centralityF};
1015          fRecMCEventHist1->Fill(vRecMCEventHist1);//
1016
1017          Double_t vRecMCEventHist2[4] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetZv()-vtxMC[2],multMBTracks, centralityF};
1018          fRecMCEventHist2->Fill(vRecMCEventHist2);
1019
1020        }//
1021
1022        //
1023        // MC histograms for track efficiency studies
1024        //
1025        for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc)
1026        {
1027          TParticle* particle = stack->Particle(iMc);
1028          if (!particle)
1029          continue;
1030
1031          Double_t vTrackMatrix[4] = {vtxMC[2],particle->Pt(),particle->Eta(),centralityF};
1032
1033          // only charged particles
1034          if(!particle->GetPDG()) continue;
1035          Double_t charge = particle->GetPDG()->Charge()/3.;
1036          if (TMath::Abs(charge) < 0.001)
1037          continue;
1038
1039          // only postive charged
1040          if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.)
1041          continue;
1042
1043          // only negative charged
1044          if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.)
1045          continue;
1046
1047          // physical primary
1048          Bool_t prim = stack->IsPhysicalPrimary(iMc);
1049
1050          // check accepted
1051          if(accCuts->AcceptTrack(particle))
1052          {
1053
1054            if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) )
1055              fGenPrimTrackMatrix->Fill(vTrackMatrix);
1056
1057            // fill control histograms
1058            if(fHistogramsOn)
1059              FillHistograms(stack,iMc,AlidNdPtHelper::kAccTracks, centralityF);
1060
1061            // check multiple found tracks
1062            Int_t multCount = 0;
1063            for(Int_t iRec=0; iRec<multRec; ++iRec)
1064            {
1065              if(iMc == labelsRec[iRec])
1066              {
1067                multCount++;
1068                if(multCount>1)
1069                {
1070                  fRecMultTrackMatrix->Fill(vTrackMatrix);
1071
1072                  // fill control histogram
1073                  if(fHistogramsOn) {
1074                    Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
1075                    Double_t vMCMultRecTrackHist1[4] = {particle->Pt(), particle->Eta(), pid, centralityF};
1076                    fMCMultRecTrackHist1->Fill(vMCMultRecTrackHist1);
1077                  }
1078                }
1079              }
1080            }
1081
1082            // check reconstructed
1083            for(Int_t iRec=0; iRec<multRec; ++iRec)
1084            {
1085              if(iMc == labelsRec[iRec])
1086              {
1087                fRecTrackMatrix->Fill(vTrackMatrix);
1088
1089                if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) {
1090                  fRecPrimTrackMatrix->Fill(vTrackMatrix);
1091                }
1092                if(!prim) fRecSecTrackMatrix->Fill(vTrackMatrix);
1093
1094                // fill control histograms
1095                if(fHistogramsOn)
1096                  FillHistograms(stack,iMc,AlidNdPtHelper::kRecTracks, centralityF);
1097
1098                break;
1099              }//if(iMc == labelsRec[iRec])
1100            }//reco tracks
1101          }//accepted tracks
1102        }//stack loop
1103      }//is triggered
1104    } // end bUseMC
1105
1106   if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;
1107   if(labelsAll) delete [] labelsAll; labelsAll = 0;
1108   if(labelsAcc) delete [] labelsAcc; labelsAcc = 0;
1109   if(labelsRec) delete [] labelsRec; labelsRec = 0;
1110
1111 }
1112
1113 //_____________________________________________________________________________
1114 void AlidNdPtAnalysisPbPb::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, Float_t centralityF) {
1115  // multiplicity  histograms
1116
1117
1118   if(!allChargedTracks) return;
1119   if(!labelsAll) return;
1120   if(!labelsAcc) return;
1121   if(!labelsRec) return;
1122
1123   Int_t entries = allChargedTracks->GetEntries();
1124   for(Int_t i=0; i<entries; ++i)
1125   {
1126      AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
1127      if(!track) continue;
1128      if(track->Charge() == 0) continue;
1129
1130      Int_t label = TMath::Abs(track->GetLabel());
1131      for(Int_t iAll=0; iAll<multAll; ++iAll) {
1132        if(label == labelsAll[iAll]) {
1133          Double_t v1[3] = {track->Pt(), multAll, centralityF};
1134          fRecTrackMultHist1[AlidNdPtHelper::kAllTracks]->Fill(v1);
1135        }
1136      }
1137      for(Int_t iAcc=0; iAcc<multAcc; ++iAcc) {
1138        if(label == labelsAcc[iAcc]) {
1139          Double_t v2[3] = {track->Pt(), multAcc, centralityF};
1140          fRecTrackMultHist1[AlidNdPtHelper::kAccTracks]->Fill(v2);
1141        }
1142      }
1143      for(Int_t iRec=0; iRec<multRec; ++iRec) {
1144        if(label == labelsRec[iRec]) {
1145          Double_t v3[3] = {track->Pt(), multRec, centralityF};
1146          fRecTrackMultHist1[AlidNdPtHelper::kRecTracks]->Fill(v3);
1147        }//out
1148      }
1149   }
1150 }
1151
1152 //_____________________________________________________________________________
1153 void AlidNdPtAnalysisPbPb::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, AlidNdPtHelper::TrackObject trackObj, Float_t centralityF)
1154 {
1155
1156   //
1157   // Fill ESD track and MC histograms
1158   //
1159   if(!esdTrack) return;
1160
1161   Float_t q = esdTrack->Charge();
1162   if(TMath::Abs(q) < 0.001) return;
1163
1164   Float_t pt = esdTrack->Pt();
1165   Float_t eta = esdTrack->Eta();
1166   Float_t phi = esdTrack->Phi();
1167
1168   Float_t dca[2], bCov[3];
1169   esdTrack->GetImpactParameters(dca,bCov);
1170
1171   Int_t nClust = esdTrack->GetTPCclusters(0);
1172   Float_t chi2PerCluster = 0.;
1173   if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);
1174
1175
1176   // fill histograms
1177   Double_t values[4] = {pt,eta,phi,centralityF};
1178   fRecTrackHist1[trackObj]->Fill(values);
1179
1180   //
1181   // Fill rec vs MC information
1182   //
1183   if(!stack) return;
1184
1185   Int_t label = TMath::Abs(esdTrack->GetLabel());
1186   //if(label == 0) return;
1187
1188   TParticle* particle = stack->Particle(label);
1189   if(!particle) return;
1190
1191   Int_t motherPdg = -1;
1192   TParticle* mother = 0;
1193
1194   //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
1195   Int_t motherLabel = particle->GetMother(0);
1196   if(motherLabel>0) mother = stack->Particle(motherLabel);
1197   if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
1198   //Int_t mech = particle->GetUniqueID(); // production mechanism
1199
1200   if(!particle->GetPDG()) return;
1201   Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
1202   if(TMath::Abs(gq)<0.001) return;
1203   Float_t gpt = particle->Pt();
1204   Float_t geta = particle->Eta();
1205
1206   Double_t dpt=0;
1207   //printf("pt %f, gpt %f \n",pt,gpt);
1208   if(gpt) dpt = (pt-gpt)/gpt;
1209   Double_t deta = (eta-geta);
1210
1211   // fill histograms
1212   if(trackObj == AlidNdPtHelper::kRecTracks)  //RecTracks???
1213   {
1214     Double_t vTrackPtCorrelationMatrix[4]={pt,gpt,geta,centralityF};
1215     fTrackPtCorrelationMatrix->Fill(vTrackPtCorrelationMatrix);
1216
1217     Double_t vRecMCTrackHist1[5]={gpt,geta,dpt,deta,centralityF};
1218     fRecMCTrackHist1->Fill(vRecMCTrackHist1);
1219   }
1220 }
1221
1222 //_____________________________________________________________________________
1223 void AlidNdPtAnalysisPbPb::FillHistograms(AliStack *const stack, Int_t label, AlidNdPtHelper::TrackObject trackObj, Float_t centralityF)
1224 {
1225
1226   // Fill MC histograms
1227   if(!stack) return;
1228
1229   TParticle* particle = stack->Particle(label);
1230   if(!particle) return;
1231
1232   Int_t motherPdg = -1;
1233   TParticle* mother = 0;
1234
1235   //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
1236   Int_t motherLabel = particle->GetMother(0);
1237   if(motherLabel>0) mother = stack->Particle(motherLabel);
1238   if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
1239   Int_t mech = particle->GetUniqueID(); // production mechanism
1240
1241   if(!particle->GetPDG()) return;
1242   Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
1243   if(TMath::Abs(gq) < 0.001) return;
1244
1245   Float_t gpt = particle->Pt();
1246   //Float_t qgpt = particle->Pt() * gq;
1247   Float_t geta = particle->Eta();
1248   Float_t gphi = particle->Phi();
1249   //Float_t gpz = particle->Pz();
1250
1251   Bool_t prim = stack->IsPhysicalPrimary(label);
1252   //Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();
1253
1254   Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);
1255
1256   //if(prim&&pid==5) printf("pdgcode %d, production mech %d \n",particle->GetPdgCode(),mech);
1257   //if(!prim) printf("motherPdg %d, particle %d, production mech %d\n",motherPdg, particle->GetPdgCode(),mech);
1258
1259   //
1260   // fill histogram
1261   //
1262   Double_t vMCTrackHist1[4] = {gpt,geta,gphi,centralityF};
1263   fMCTrackHist1[trackObj]->Fill(vMCTrackHist1);
1264
1265   Double_t vMCPrimTrackHist1[6] = {gpt,geta,pid,mech,motherPdg,centralityF};
1266   Double_t vMCPrimTrackHist2[4] = {TMath::Abs(particle->GetPdgCode()),mech,motherPdg,centralityF};
1267
1268   //if(prim && AliPWG0Helper::IsPrimaryCharged(particle, label)) fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
1269
1270   if(prim) {
1271     fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
1272     if(pid == 5) fMCPrimTrackHist2[trackObj]->Fill(vMCPrimTrackHist2);
1273   }
1274   else {
1275     fMCSecTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);
1276   }
1277
1278 }
1279
1280 //_____________________________________________________________________________
1281 Long64_t AlidNdPtAnalysisPbPb::Merge(TCollection* const list)
1282 {
1283   // Merge list of objects (needed by PROOF)
1284
1285   //  init if not done already
1286   if (!fIsInit) { Init(); }
1287
1288   if (!list)
1289   return 0;
1290
1291   if (list->IsEmpty())
1292   return 1;
1293
1294   TIterator* iter = list->MakeIterator();
1295   TObject* obj = 0;
1296
1297   //
1298   //TList *collPhysSelection = new TList;
1299
1300   // collection of generated histograms
1301
1302   Int_t count=0;
1303   while((obj = iter->Next()) != 0) {
1304     AlidNdPtAnalysisPbPb* entry = dynamic_cast<AlidNdPtAnalysisPbPb*>(obj);
1305     if (entry == 0) continue;
1306
1307     // physics selection
1308     //printf("entry->GetPhysicsTriggerSelection() %p \n", entry->GetPhysicsTriggerSelection());
1309     //AliPhysicsSelection *physSel = entry->GetPhysicsTriggerSelection();
1310     //if( physSel ) collPhysSelection->Add(physSel);
1311
1312     //
1313     fTrackPtCorrelationMatrix->Add(entry->fTrackPtCorrelationMatrix);
1314
1315     //
1316     fGenEventMatrix->Add(entry->fGenEventMatrix);
1317
1318     fTriggerEventMatrix->Add(entry->fTriggerEventMatrix);
1319
1320     fRecEventMatrix->Add(entry->fRecEventMatrix);
1321     //
1322     fGenTrackEventMatrix->Add(entry->fGenTrackEventMatrix);
1323
1324     fTriggerTrackEventMatrix->Add(entry->fTriggerTrackEventMatrix);
1325
1326     fRecTrackEventMatrix->Add(entry->fRecTrackEventMatrix);
1327
1328     //
1329     fGenTrackMatrix->Add(entry->fGenTrackMatrix);
1330     fGenPrimTrackMatrix->Add(entry->fGenPrimTrackMatrix);
1331     fRecPrimTrackMatrix->Add(entry->fRecPrimTrackMatrix);
1332     //
1333     fRecTrackMatrix->Add(entry->fRecTrackMatrix);
1334     fRecSecTrackMatrix->Add(entry->fRecSecTrackMatrix);
1335     //
1336     fRecMultTrackMatrix->Add(entry->fRecMultTrackMatrix);
1337
1338     //
1339     // control analysis histograms
1340     //
1341     fMCEventHist1->Add(entry->fMCEventHist1);
1342     fRecEventHist1->Add(entry->fRecEventHist1);
1343     fRecEventHist2->Add(entry->fRecEventHist2);
1344     fRecMCEventHist1->Add(entry->fRecMCEventHist1);
1345     fRecMCEventHist2->Add(entry->fRecMCEventHist2);
1346
1347
1348     for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
1349       fMCTrackHist1[i]->Add(entry->fMCTrackHist1[i]);
1350
1351       fMCPrimTrackHist1[i]->Add(entry->fMCPrimTrackHist1[i]);
1352       fMCPrimTrackHist2[i]->Add(entry->fMCPrimTrackHist2[i]);
1353       fMCSecTrackHist1[i]->Add(entry->fMCSecTrackHist1[i]);
1354
1355       fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);
1356       fRecTrackHist2[i]->Add(entry->fRecTrackHist2[i]);
1357       fRecTrackMultHist1[i]->Add(entry->fRecTrackMultHist1[i]);
1358     }
1359     fRecMCTrackHist1->Add(entry->fRecMCTrackHist1);
1360     fMCMultRecTrackHist1->Add(entry->fMCMultRecTrackHist1);
1361     fRecTrackHist3->Add(entry->fRecTrackHist3);
1362
1363   count++;
1364   }
1365
1366   //AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();
1367   //if( trigSelection ) trigSelection->Merge(collPhysSelection);
1368   //if(collPhysSelection) delete collPhysSelection;
1369
1370 return count;
1371 }
1372
1373
1374
1375 //_____________________________________________________________________________
1376 void AlidNdPtAnalysisPbPb::Analyse()
1377 {
1378
1379   //  init if not done already
1380   if (!fIsInit) { Init(); }
1381
1382   // Analyse histograms
1383   //
1384   TH1::AddDirectory(kFALSE);
1385   TH1 *h=0, *h1=0, *h2=0, *h2c = 0;
1386   THnSparse *hs=0;
1387   TH2 *h2D=0;
1388
1389   char name[256];
1390   TObjArray *aFolderObj = new TObjArray;
1391
1392   //
1393   // LHC backgraund in all and 0-bins
1394   //
1395   //AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection();
1396   //trigSel->SaveHistograms("physics_selection");
1397
1398   //
1399   // Reconstructed event vertex
1400   //
1401   h = fRecEventHist1->Projection(0);
1402   h->SetName("Xv");
1403   aFolderObj->Add(h);
1404
1405   h = fRecEventHist1->Projection(1);
1406   h->SetName("Yv");
1407   aFolderObj->Add(h);
1408
1409   h = fRecEventHist1->Projection(2);
1410   h->SetName("Zv");
1411   aFolderObj->Add(h);
1412
1413   //
1414   // multiplicity
1415   //
1416   h = fRecEventHist2->Projection(1);
1417   h->SetName("multMB");
1418   aFolderObj->Add(h);
1419
1420   h2D = fRecEventHist2->Projection(0,1);
1421   h2D->SetName("Zv_vs_multiplicity_MB");
1422   aFolderObj->Add(h2D);
1423
1424   //
1425   // reconstructed pt histograms
1426   //
1427   h = fRecTrackHist1[0]->Projection(0);
1428   h->Scale(1.,"width");
1429   h->SetName("pt_all_ch");
1430   aFolderObj->Add(h);
1431
1432   h = fRecTrackHist1[1]->Projection(0);
1433   h->Scale(1.,"width");
1434   h->SetName("pt_acc");
1435   aFolderObj->Add(h);
1436
1437   h = fRecTrackHist1[2]->Projection(0);
1438   h->Scale(1.,"width");
1439   h->SetName("pt_rec");
1440   aFolderObj->Add(h);
1441
1442   //
1443   // reconstructed eta histograms
1444   //
1445   h = fRecTrackHist1[0]->Projection(1);
1446   h->SetName("eta_all_ch");
1447   aFolderObj->Add(h);
1448
1449   h = fRecTrackHist1[1]->Projection(1);
1450   h->SetName("eta_acc");
1451   aFolderObj->Add(h);
1452
1453   h = fRecTrackHist1[2]->Projection(1);
1454   h->SetName("eta_rec");
1455   aFolderObj->Add(h);
1456
1457   //
1458   // reconstructed phi histograms
1459   //
1460   h = fRecTrackHist1[0]->Projection(2);
1461   h->SetName("phi_all_ch");
1462   aFolderObj->Add(h);
1463
1464   h = fRecTrackHist1[1]->Projection(2);
1465   h->SetName("phi_acc");
1466   aFolderObj->Add(h);
1467
1468   h = fRecTrackHist1[2]->Projection(2);
1469   h->SetName("phi_rec");
1470   aFolderObj->Add(h);
1471
1472   //
1473   // reconstructed eta:pt histograms
1474   //
1475   h2D = fRecTrackHist1[0]->Projection(1,0);
1476   h2D->SetName("pt_eta_all_ch");
1477   aFolderObj->Add(h2D);
1478
1479   h2D = fRecTrackHist1[1]->Projection(1,0);
1480   h2D->SetName("pt_eta_acc");
1481   aFolderObj->Add(h2D);
1482
1483   h2D = fRecTrackHist1[2]->Projection(1,0);
1484   h2D->SetName("pt_eta_rec");
1485   aFolderObj->Add(h2D);
1486
1487   //
1488   // reconstructed phi:pt histograms
1489   //
1490   h2D = fRecTrackHist1[0]->Projection(2,0);
1491   h2D->SetName("pt_phi_all_ch");
1492   aFolderObj->Add(h2D);
1493
1494   h2D = fRecTrackHist1[1]->Projection(2,0);
1495   h2D->SetName("pt_phi_acc");
1496   aFolderObj->Add(h2D);
1497
1498   h2D = fRecTrackHist1[2]->Projection(2,0);
1499   h2D->SetName("pt_phi_rec");
1500   aFolderObj->Add(h2D);
1501
1502   //
1503   // reconstructed phi:eta histograms
1504   //
1505   h2D = fRecTrackHist1[0]->Projection(2,1);
1506   h2D->SetName("eta_phi_all_ch");
1507   aFolderObj->Add(h2D);
1508
1509   h2D = fRecTrackHist1[1]->Projection(2,1);
1510   h2D->SetName("eta_phi_acc");
1511   aFolderObj->Add(h2D);
1512
1513   h2D = fRecTrackHist1[2]->Projection(2,1);
1514   h2D->SetName("eta_phi_rec");
1515   aFolderObj->Add(h2D);
1516
1517   //
1518   // reconstructed nClust, chi2 vs pt, eta, phi
1519   //
1520   if(fHistogramsOn) {
1521
1522     h2D = fRecTrackHist3->Projection(0,1);
1523     h2D->SetName("nClust_chi2_rec");
1524     aFolderObj->Add(h2D);
1525
1526     h2D = fRecTrackHist3->Projection(0,2);
1527     h2D->SetName("nClust_pt_rec");
1528     aFolderObj->Add(h2D);
1529
1530     h2D = fRecTrackHist3->Projection(0,3);
1531     h2D->SetName("nClust_eta_rec");
1532     aFolderObj->Add(h2D);
1533
1534     h2D = fRecTrackHist3->Projection(0,4);
1535     h2D->SetName("nClust_phi_rec");
1536     aFolderObj->Add(h2D);
1537
1538     h2D = fRecTrackHist3->Projection(1,2);
1539     h2D->SetName("chi2_pt_rec");
1540     aFolderObj->Add(h2D);
1541
1542     h2D = fRecTrackHist3->Projection(1,3);
1543     h2D->SetName("chi2_eta_rec");
1544     aFolderObj->Add(h2D);
1545
1546     h2D = fRecTrackHist3->Projection(1,4);
1547     h2D->SetName("chi2_phi_rec");
1548     aFolderObj->Add(h2D);
1549
1550   }
1551
1552   //
1553   // calculate corrections for empty events
1554   // with multMB==0
1555   //
1556
1557   //
1558   // normalised zv to generate zv for triggered events
1559   //
1560   h = fRecEventHist2->Projection(0);
1561   if( h->Integral() ) h->Scale(1./h->Integral());
1562   h->SetName("zv_distribution_norm");
1563   aFolderObj->Add(h);
1564
1565   //
1566   // MC available
1567   //
1568   if(IsUseMCInfo()) {
1569
1570   //
1571   // Event vertex resolution
1572   //
1573   h2D = fRecMCEventHist2->Projection(0,2);
1574   h2D->SetName("DeltaXv_vs_mult");
1575   aFolderObj->Add(h2D);
1576
1577   h2D = fRecMCEventHist2->Projection(1,2);
1578   h2D->SetName("DeltaZv_vs_mult");
1579   aFolderObj->Add(h2D);
1580
1581   //
1582   // normalised zv to get trigger/trigger+vertex event differences
1583   // F(zv) = E_trig(zv,0)/Int(E_trig(zv,0) / Sum(E_trigvtx(zv,n))/Sum(Int(E_trigvtx(zv,n))dzv)
1584   //
1585   fTriggerEventMatrix->GetAxis(1)->SetRangeUser(0.,0.);
1586   h = fTriggerEventMatrix->Projection(0);
1587   h2D = fTriggerEventMatrix->Projection(0,1);
1588   if(h2D->Integral()) h->Scale(1./h2D->Integral());
1589
1590   h1 = fRecEventMatrix->Projection(0);
1591   h2D = fRecEventMatrix->Projection(0,1);
1592   if(h2D->Integral()) h1->Scale(1./h2D->Integral());
1593
1594   h->Divide(h1);
1595   h->SetName("zv_empty_events_norm");
1596   aFolderObj->Add(h);
1597
1598   fTriggerEventMatrix->GetAxis(1)->SetRange(1,fTriggerEventMatrix->GetAxis(1)->GetNbins());
1599
1600   //
1601   // rec. vs true track pt correlation matrix
1602   //
1603   hs = (THnSparse*)fTrackPtCorrelationMatrix->Clone("track_pt_correlation_matrix");
1604   aFolderObj->Add(hs);
1605
1606   //
1607   // trigger efficiency for INEL
1608   //
1609   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fGenEventMatrix->Projection(0),"zv_trig_INEL_eff_matrix");
1610   aFolderObj->Add(h);
1611
1612
1613   //
1614   // trigger bias correction (MB to INEL)
1615   //
1616   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoInel_corr_matrix");
1617   aFolderObj->Add(hs);
1618
1619   h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoInel_corr_matrix");
1620   aFolderObj->Add(h);
1621
1622   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0,1),fTriggerEventMatrix->Projection(0,1),"zv_mult_trig_MBtoInel_corr_matrix_2D");
1623   aFolderObj->Add(h2D);
1624
1625
1626   h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoInel_corr_matrix");
1627   aFolderObj->Add(h);
1628
1629
1630   //
1631   // event vertex reconstruction correction (MB)
1632   //
1633   hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix,fRecEventMatrix,"zv_mult_event_corr_matrix");
1634   aFolderObj->Add(hs);
1635
1636   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0,1),fRecEventMatrix->Projection(0,1),"zv_mult_event_corr_matrix_2D");
1637   aFolderObj->Add(h2D);
1638
1639
1640   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(1),fRecEventMatrix->Projection(1),"mult_event_corr_matrix");
1641   aFolderObj->Add(h);
1642
1643
1644   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fRecEventMatrix->Projection(0),"zv_event_corr_matrix");
1645   aFolderObj->Add(h);
1646
1647
1648
1649   //
1650   // track-event trigger bias correction (MB to INEL)
1651   //
1652   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoInel_corr_matrix");
1653   aFolderObj->Add(hs);
1654
1655   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoInel_corr_matrix");
1656   aFolderObj->Add(h2D);
1657
1658   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoInel_corr_matrix");
1659   aFolderObj->Add(h2D);
1660
1661   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoInel_corr_matrix");
1662   aFolderObj->Add(h2D);
1663
1664   // efficiency
1665
1666   h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1),fGenTrackEventMatrix->Projection(1),"pt_track_trig_MBtoInel_eff_matrix");
1667   aFolderObj->Add(h);
1668
1669
1670   //
1671   // track-event vertex reconstruction correction (MB)
1672   //
1673   hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix,fRecTrackEventMatrix,"zv_pt_eta_track_event_corr_matrix");
1674   aFolderObj->Add(hs);
1675
1676   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,2),fRecTrackEventMatrix->Projection(1,2),"eta_pt_track_event_corr_matrix");
1677   aFolderObj->Add(h2D);
1678
1679   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,0),fRecTrackEventMatrix->Projection(1,0),"pt_zv_track_event_corr_matrix");
1680   aFolderObj->Add(h2D);
1681
1682   h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(2,0),fRecTrackEventMatrix->Projection(2,0),"zv_eta_track_event_corr_matrix");
1683   aFolderObj->Add(h2D);
1684
1685   // efficiency
1686
1687   h = AlidNdPtHelper::GenerateCorrMatrix(fRecTrackEventMatrix->Projection(1),fTriggerTrackEventMatrix->Projection(1),"pt_track_event_eff_matrix");
1688   aFolderObj->Add(h);
1689
1690
1691   //
1692   // track rec. efficiency correction
1693   //
1694   hs = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix,fRecPrimTrackMatrix,"zv_pt_eta_track_corr_matrix");
1695   aFolderObj->Add(hs);
1696
1697   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,2),fRecPrimTrackMatrix->Projection(1,2),"eta_pt_track_corr_matrix");
1698   aFolderObj->Add(h2D);
1699
1700   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,0),fRecPrimTrackMatrix->Projection(1,0),"pt_zv_track_corr_matrix");
1701   aFolderObj->Add(h2D);
1702
1703   h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2,0),fRecPrimTrackMatrix->Projection(2,0),"zv_eta_track_corr_matrix");
1704   aFolderObj->Add(h2D);
1705
1706
1707   h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(0),fRecPrimTrackMatrix->Projection(0),"zv_track_corr_matrix");
1708   aFolderObj->Add(h);
1709
1710   h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1),fRecPrimTrackMatrix->Projection(1),"pt_track_corr_matrix");
1711   aFolderObj->Add(h);
1712
1713   // efficiency
1714
1715   h = AlidNdPtHelper::GenerateCorrMatrix(fRecPrimTrackMatrix->Projection(1), fGenPrimTrackMatrix->Projection(1),"pt_track_eff_matrix");
1716   aFolderObj->Add(h);
1717
1718   h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2),fRecPrimTrackMatrix->Projection(2),"eta_track_corr_matrix");
1719   aFolderObj->Add(h);
1720
1721   //
1722   // secondary track contamination correction
1723   //
1724   //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");
1725   hs = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");
1726   aFolderObj->Add(hs);
1727
1728   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_track_cont_matrix");
1729   aFolderObj->Add(h2D);
1730
1731   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_track_cont_matrix");
1732   aFolderObj->Add(h2D);
1733
1734   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_track_cont_matrix");
1735   aFolderObj->Add(h2D);
1736
1737   h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_track_cont_matrix");
1738   aFolderObj->Add(h);
1739
1740
1741   h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_track_cont_matrix");
1742   aFolderObj->Add(h);
1743
1744   h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_track_cont_matrix");
1745   aFolderObj->Add(h);
1746
1747   //
1748   // multiple track reconstruction correction
1749   //
1750   //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");
1751   hs = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");
1752   aFolderObj->Add(hs);
1753
1754   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_mult_track_cont_matrix");
1755   aFolderObj->Add(h2D);
1756
1757   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_mult_track_cont_matrix");
1758   aFolderObj->Add(h2D);
1759
1760   h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_mult_track_cont_matrix");
1761   aFolderObj->Add(h2D);
1762
1763   h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_mult_track_cont_matrix");
1764   aFolderObj->Add(h);
1765
1766   h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_mult_track_cont_matrix");
1767   aFolderObj->Add(h);
1768
1769   h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_mult_track_cont_matrix");
1770   aFolderObj->Add(h);
1771
1772   //
1773   // Control histograms
1774   //
1775
1776   if(fHistogramsOn) {
1777
1778   // Efficiency electrons, muons, pions, kaons, protons, all
1779   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,1);
1780   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,1);
1781   h1 = fMCPrimTrackHist1[1]->Projection(0);
1782   h2 = fMCPrimTrackHist1[2]->Projection(0);
1783   h2c = (TH1D *)h2->Clone();
1784   h2c->Divide(h1);
1785   h2c->SetName("eff_pt_electrons");
1786   aFolderObj->Add(h2c);
1787
1788   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(2,2);
1789   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(2,2);
1790   h1 = fMCPrimTrackHist1[1]->Projection(0);
1791   h2 = fMCPrimTrackHist1[2]->Projection(0);
1792   h2c = (TH1D *)h2->Clone();
1793   h2c->Divide(h1);
1794   h2c->SetName("eff_pt_muons");
1795   aFolderObj->Add(h2c);
1796
1797   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(3,3);
1798   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(3,3);
1799   h1 = fMCPrimTrackHist1[1]->Projection(0);
1800   h2 = fMCPrimTrackHist1[2]->Projection(0);
1801   h2c = (TH1D *)h2->Clone();
1802   h2c->Divide(h1);
1803   h2c->SetName("eff_pt_pions");
1804   aFolderObj->Add(h2c);
1805
1806   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(4,4);
1807   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(4,4);
1808   h1 = fMCPrimTrackHist1[1]->Projection(0);
1809   h2 = fMCPrimTrackHist1[2]->Projection(0);
1810   h2c = (TH1D *)h2->Clone();
1811   h2c->Divide(h1);
1812   h2c->SetName("eff_pt_kaons");
1813   aFolderObj->Add(h2c);
1814
1815   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(5,5);
1816   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(5,5);
1817   h1 = fMCPrimTrackHist1[1]->Projection(0);
1818   h2 = fMCPrimTrackHist1[2]->Projection(0);
1819   h2c = (TH1D *)h2->Clone();
1820   h2c->Divide(h1);
1821   h2c->SetName("eff_pt_protons");
1822   aFolderObj->Add(h2c);
1823
1824   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,5);
1825   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,5);
1826   h1 = fMCPrimTrackHist1[1]->Projection(0);
1827   h2 = fMCPrimTrackHist1[2]->Projection(0);
1828   h2c = (TH1D *)h2->Clone();
1829   h2c->Divide(h1);
1830   h2c->SetName("eff_pt_selected");
1831   aFolderObj->Add(h2c);
1832
1833   fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,6);
1834   fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,6);
1835   h1 = fMCPrimTrackHist1[1]->Projection(0);
1836   h2 = fMCPrimTrackHist1[2]->Projection(0);
1837   h2c = (TH1D *)h2->Clone();
1838   h2c->Divide(h1);
1839   h2c->SetName("eff_pt_all");
1840   aFolderObj->Add(h2c);
1841
1842   fMCPrimTrackHist1[1]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[1]->GetAxis(1)->GetNbins());
1843   fMCPrimTrackHist1[2]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[2]->GetAxis(1)->GetNbins());
1844
1845   //  pt spetra
1846   // - rec, primaries, secondaries
1847   // - primaries (pid)
1848   // - secondaries (pid)
1849   // - secondaries (mech)
1850   // - secondaries (mother)
1851   //
1852
1853   TH1D *mcPtAccall = fMCTrackHist1[1]->Projection(0);
1854   mcPtAccall->SetName("mc_pt_acc_all");
1855   aFolderObj->Add(mcPtAccall);
1856
1857   TH1D *mcPtAccprim = fMCPrimTrackHist1[1]->Projection(0);
1858   mcPtAccprim->SetName("mc_pt_acc_prim");
1859   aFolderObj->Add(mcPtAccprim);
1860
1861   TH1D *mcPtRecall = fMCTrackHist1[2]->Projection(0);
1862   mcPtRecall->SetName("mc_pt_rec_all");
1863   aFolderObj->Add(mcPtRecall);
1864
1865   TH1D *mcPtRecprim = fMCPrimTrackHist1[2]->Projection(0);
1866   mcPtRecprim->SetName("mc_pt_rec_prim");
1867   aFolderObj->Add(mcPtRecprim);
1868
1869   TH1D *mcPtRecsec = fMCSecTrackHist1[2]->Projection(0);
1870   mcPtRecsec->SetName("mc_pt_rec_sec");
1871   aFolderObj->Add(mcPtRecsec);
1872
1873   for(Int_t i = 0; i<6; i++)
1874   {
1875     sprintf(name,"mc_pt_rec_prim_pid_%d",i);
1876     fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
1877     h = fMCPrimTrackHist1[2]->Projection(0);
1878     h->SetName(name);
1879     aFolderObj->Add(h);
1880
1881     sprintf(name,"mc_pt_rec_sec_pid_%d",i);
1882     fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
1883     h = fMCSecTrackHist1[2]->Projection(0);
1884     h->SetName(name);
1885     aFolderObj->Add(h);
1886
1887     // production mechanisms for given pid
1888     fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);
1889
1890     for(Int_t j=0; j<20; j++) {
1891       if(j == 4) {
1892         // decay
1893
1894         sprintf(name,"mc_pt_rec_sec_pid_%d_decay",i);
1895         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1896         h = fMCSecTrackHist1[2]->Projection(0);
1897         h->SetName(name);
1898         aFolderObj->Add(h);
1899
1900         sprintf(name,"mc_eta_rec_sec_pid_%d_decay",i);
1901         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1902         h = fMCSecTrackHist1[2]->Projection(1);
1903         h->SetName(name);
1904         aFolderObj->Add(h);
1905
1906         sprintf(name,"mc_mother_rec_sec_pid_%d_decay",i);
1907         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1908         h = fMCSecTrackHist1[2]->Projection(4);
1909         h->SetName(name);
1910         aFolderObj->Add(h);
1911
1912       } else if (j == 5) {
1913        // conversion
1914
1915         sprintf(name,"mc_pt_rec_sec_pid_%d_conv",i);
1916         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1917         h = fMCSecTrackHist1[2]->Projection(0);
1918         h->SetName(name);
1919         aFolderObj->Add(h);
1920
1921         sprintf(name,"mc_eta_rec_sec_pid_%d_conv",i);
1922         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1923         h = fMCSecTrackHist1[2]->Projection(1);
1924         h->SetName(name);
1925         aFolderObj->Add(h);
1926
1927         sprintf(name,"mc_mother_rec_sec_pid_%d_conv",i);
1928         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1929         h = fMCSecTrackHist1[2]->Projection(4);
1930         h->SetName(name);
1931         aFolderObj->Add(h);
1932
1933      } else if (j == 13) {
1934        // mat
1935
1936         sprintf(name,"mc_pt_rec_sec_pid_%d_mat",i);
1937         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1938         h = fMCSecTrackHist1[2]->Projection(0);
1939         h->SetName(name);
1940         aFolderObj->Add(h);
1941
1942         sprintf(name,"mc_eta_rec_sec_pid_%d_mat",i);
1943         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1944         h = fMCSecTrackHist1[2]->Projection(1);
1945         h->SetName(name);
1946         aFolderObj->Add(h);
1947
1948         sprintf(name,"mc_eta_mother_rec_sec_pid_%d_mat",i);
1949         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1950         h = fMCSecTrackHist1[2]->Projection(4,1);
1951         h->SetName(name);
1952         aFolderObj->Add(h);
1953
1954         sprintf(name,"mc_mother_rec_sec_pid_%d_mat",i);
1955         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1956         h = fMCSecTrackHist1[2]->Projection(4);
1957         h->SetName(name);
1958         aFolderObj->Add(h);
1959
1960         sprintf(name,"mc_pt_mother_rec_sec_pid_%d_mat",i);
1961         fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);
1962         h = fMCSecTrackHist1[2]->Projection(4,0);
1963         h->SetName(name);
1964         aFolderObj->Add(h);
1965
1966      } else {
1967        continue;
1968      }
1969    }
1970
1971   }
1972   } // end fHistogramOn
1973
1974   //
1975   //  resolution histograms
1976   //  only for reconstructed tracks
1977   //
1978
1979   TH2F *h2F=0;
1980   TCanvas * c = new TCanvas("resol","resol");
1981   c->cd();
1982
1983   //
1984   fRecMCTrackHist1->GetAxis(1)->SetRangeUser(-0.8,0.79);
1985
1986   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);
1987   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
1988   h->SetXTitle("p_{tmc} (GeV/c)");
1989   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");
1990   h->Draw();
1991   h->SetName("pt_resolution_vs_mcpt");
1992   aFolderObj->Add(h);
1993
1994   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);
1995   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
1996   h->SetXTitle("p_{tmc} (GeV/c)");
1997   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");
1998   h->Draw();
1999   h->SetName("dpt_mean_vs_mcpt");
2000   aFolderObj->Add(h);
2001
2002   //
2003   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);
2004   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
2005   h->SetXTitle("p_{tmc} (GeV/c)");
2006   h->SetYTitle("(#eta-#eta_{mc}) resolution");
2007   h->Draw();
2008   h->SetName("eta_resolution_vs_mcpt");
2009   aFolderObj->Add(h);
2010
2011   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);
2012   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
2013   h->SetXTitle("p_{tmc} (GeV/c)");
2014   h->SetYTitle("(#eta-mc#eta) mean");
2015   h->Draw();
2016   h->SetName("deta_mean_vs_mcpt");
2017   aFolderObj->Add(h);
2018
2019   //
2020   fRecMCTrackHist1->GetAxis(1)->SetRange(1,fRecMCTrackHist1->GetAxis(1)->GetNbins());
2021
2022   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);
2023   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
2024   h->SetXTitle("#eta_{mc}");
2025   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");
2026   h->Draw();
2027   h->SetName("pt_resolution_vs_mceta");
2028   aFolderObj->Add(h);
2029
2030   h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);
2031   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
2032   h->SetXTitle("#eta_{mc}");
2033   h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");
2034   h->Draw();
2035   h->SetName("dpt_mean_vs_mceta");
2036   aFolderObj->Add(h);
2037
2038   //
2039   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);
2040   h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);
2041   h->SetXTitle("#eta_{mc}");
2042   h->SetYTitle("(#eta-#eta_{mc}) resolution");
2043   h->Draw();
2044   h->SetName("eta_resolution_vs_mceta");
2045   aFolderObj->Add(h);
2046
2047   h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);
2048   h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);
2049   h->SetXTitle("#eta_{mc}");
2050   h->SetYTitle("(#eta-mc#eta) mean");
2051   h->Draw();
2052   h->SetName("deta_mean_vs_mceta");
2053   aFolderObj->Add(h);
2054
2055   fRecMCTrackHist1->GetAxis(0)->SetRange(1,fRecMCTrackHist1->GetAxis(0)->GetNbins());
2056
2057   } // end use MC info
2058
2059   // export objects to analysis folder
2060   fAnalysisFolder = ExportToFolder(aFolderObj);
2061
2062   // delete only TObjArray
2063   if(aFolderObj) delete aFolderObj;
2064 }
2065
2066 //_____________________________________________________________________________
2067 TFolder* AlidNdPtAnalysisPbPb::ExportToFolder(TObjArray * const array)
2068 {
2069   // recreate folder avery time and export objects to new one
2070   //
2071   AlidNdPtAnalysisPbPb * comp=this;
2072   TFolder *folder = comp->GetAnalysisFolder();
2073
2074   TString name, title;
2075   TFolder *newFolder = 0;
2076   Int_t i = 0;
2077   Int_t size = array->GetSize();
2078
2079   if(folder) {
2080      // get name and title from old folder
2081      name = folder->GetName();
2082      title = folder->GetTitle();
2083
2084          // delete old one
2085      delete folder;
2086
2087          // create new one
2088      newFolder = CreateFolder(name.Data(),title.Data());
2089      newFolder->SetOwner();
2090
2091          // add objects to folder
2092      while(i < size) {
2093            newFolder->Add(array->At(i));
2094            i++;
2095          }
2096   }
2097
2098 return newFolder;
2099 }
2100
2101 //_____________________________________________________________________________
2102 TFolder* AlidNdPtAnalysisPbPb::CreateFolder(TString name,TString title) {
2103 // create folder for analysed histograms
2104 //
2105 TFolder *folder = 0;
2106   folder = new TFolder(name.Data(),title.Data());
2107
2108   return folder;
2109 }