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