pT analysis in PbPb collisions (Markus Koehler)
authorjotwinow <jotwinow@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 21 Oct 2010 17:05:51 +0000 (17:05 +0000)
committerjotwinow <jotwinow@f7af4fe6-9843-0410-8265-dc069ae4e863>
Thu, 21 Oct 2010 17:05:51 +0000 (17:05 +0000)
PWG0/PWG0depLinkDef.h
PWG0/dNdPt/AlidNdPtAnalysisPbPb.cxx [new file with mode: 0644]
PWG0/dNdPt/AlidNdPtAnalysisPbPb.h [new file with mode: 0644]
PWG0/libPWG0dep.pkg

index bc132cd..b4f2190 100644 (file)
@@ -10,6 +10,7 @@
 
 #pragma link C++ class AlidNdPtTask+;
 #pragma link C++ class AlidNdPtHelper+;
+#pragma link C++ class AlidNdPtAnalysisPbPb+;
 #pragma link C++ class AlidNdPtAnalysis+;
 #pragma link C++ class AlidNdPtCorrection+;
 
diff --git a/PWG0/dNdPt/AlidNdPtAnalysisPbPb.cxx b/PWG0/dNdPt/AlidNdPtAnalysisPbPb.cxx
new file mode 100644 (file)
index 0000000..744b5fa
--- /dev/null
@@ -0,0 +1,1949 @@
+/**************************************************************************\r
+ * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
+ *                                                                        *\r
+ * Author: The ALICE Off-line Project.                                    *\r
+ * Contributors are mentioned in the code where appropriate.              *\r
+ *                                                                        *\r
+ * Permission to use, copy, modify and distribute this software and its   *\r
+ * documentation strictly for non-commercial purposes is hereby granted   *\r
+ * without fee, provided that the above copyright notice appears in all   *\r
+ * copies and that both the copyright notice and this permission notice   *\r
+ * appear in the supporting documentation. The authors make no claims     *\r
+ * about the suitability of this software for any purpose. It is          *\r
+ * provided "as is" without express or implied warranty.                  *\r
+ **************************************************************************/\r
+//------------------------------------------------------------------------------\r
+// AlidNdPtAnalysisPbPb class. \r
+// \r
+// a. functionality:\r
+// - fills analysis control histograms\r
+// - fills generic correction matrices \r
+// - generates correction matrices \r
+//\r
+// b. data members:\r
+// - generic correction matrices\r
+// - control histograms\r
+//\r
+// Author: J.Otwinowski 04/11/2008 \r
+//------------------------------------------------------------------------------\r
+\r
+#include "TH1.h"\r
+#include "TH2.h"\r
+#include "TCanvas.h"\r
+#include "THnSparse.h"\r
+\r
+#include "AliHeader.h"  \r
+#include "AliGenEventHeader.h"  \r
+#include "AliStack.h"  \r
+#include "AliESDEvent.h"  \r
+#include "AliMCEvent.h"  \r
+#include "AliESDtrackCuts.h"  \r
+#include "AliLog.h" \r
+#include "AliMultiplicity.h"\r
+#include "AliTracker.h"\r
+\r
+#include "AlidNdPtEventCuts.h"\r
+#include "AlidNdPtAcceptanceCuts.h"\r
+#include "AliPhysicsSelection.h"\r
+#include "AliTriggerAnalysis.h"\r
+\r
+#include "AliPWG0Helper.h"\r
+#include "AlidNdPtHelper.h"\r
+#include "AlidNdPtAnalysisPbPb.h"\r
+\r
+\r
+using namespace std;\r
+\r
+ClassImp(AlidNdPtAnalysisPbPb)\r
+\r
+//_____________________________________________________________________________\r
+  AlidNdPtAnalysisPbPb::AlidNdPtAnalysisPbPb(): AlidNdPt(),\r
+  fAnalysisFolder(0),\r
+  fHistogramsOn(kFALSE),\r
+\r
+  // rec. track pt vs true track pt correlation matrix \r
+  fTrackPtCorrelationMatrix(0),\r
+\r
+  // event level correction\r
+  fGenEventMatrix(0),\r
+  fTriggerEventMatrix(0),\r
+  fRecEventMatrix(0),\r
+\r
+  //\r
+  // track-event level correction \r
+  //\r
+  fGenTrackEventMatrix(0),\r
+\r
+  fTriggerTrackEventMatrix(0),\r
+\r
+  fRecTrackEventMatrix(0),\r
+\r
+  // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)\r
+  fGenTrackMatrix(0),\r
+  fGenPrimTrackMatrix(0),\r
+  fRecPrimTrackMatrix(0),\r
+\r
+  // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)\r
+  fRecTrackMatrix(0),\r
+  fRecSecTrackMatrix(0),\r
+\r
+  // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix)\r
+  fRecMultTrackMatrix(0),\r
+\r
+  // event control histograms\r
+  fMCEventHist1(0),\r
+  fRecEventHist1(0),\r
+  fRecEventHist2(0),\r
+  fRecMCEventHist1(0),\r
+  fRecMCEventHist2(0),\r
+\r
+  // rec. pt and eta resolution w.r.t MC\r
+  fRecMCTrackHist1(0),\r
+\r
+  //multple reconstructed tracks\r
+  fMCMultRecTrackHist1(0), \r
+\r
+  // rec. track control histograms\r
+  fRecTrackHist2(0)\r
+{\r
+  // default constructor\r
+  for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { \r
+    fMCTrackHist1[i]=0;     \r
+    fMCPrimTrackHist1[i]=0;     \r
+    fMCPrimTrackHist2[i]=0;     \r
+    fMCSecTrackHist1[i]=0;     \r
+    fRecTrackHist1[i]=0;     \r
+    fRecTrackMultHist1[i]=0;     \r
+  }\r
+  Init();\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AlidNdPtAnalysisPbPb::AlidNdPtAnalysisPbPb(Char_t* name, Char_t* title): AlidNdPt(name,title),\r
+  fAnalysisFolder(0),\r
+  fHistogramsOn(kFALSE),\r
+\r
+  // rec. track pt vs true track pt correlation matrix \r
+  fTrackPtCorrelationMatrix(0),\r
+\r
+  // event level correction\r
+  fGenEventMatrix(0),\r
+\r
+  fTriggerEventMatrix(0),\r
+\r
+  fRecEventMatrix(0),\r
+\r
+  //\r
+  // track-event level correction \r
+  //\r
+  fGenTrackEventMatrix(0),\r
+\r
+  fTriggerTrackEventMatrix(0),\r
+\r
+  fRecTrackEventMatrix(0),\r
+\r
+  // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)\r
+  fGenTrackMatrix(0),\r
+  fGenPrimTrackMatrix(0),\r
+  fRecPrimTrackMatrix(0),\r
+\r
+  // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)\r
+  fRecTrackMatrix(0),\r
+  fRecSecTrackMatrix(0),\r
+\r
+  // multiple rec. track contamination corrections (fRecMultTrackMatrix / fRecTrackMatrix)\r
+  fRecMultTrackMatrix(0),\r
+\r
+  // event control histograms\r
+  fMCEventHist1(0),\r
+  fRecEventHist1(0),\r
+  fRecEventHist2(0),\r
+  fRecMCEventHist1(0),\r
+  fRecMCEventHist2(0),\r
\r
+  // rec. pt and eta resolution w.r.t MC\r
+  fRecMCTrackHist1(0),\r
+\r
+  //multple reconstructed tracks\r
+  fMCMultRecTrackHist1(0), \r
+\r
+  // rec. track control histograms\r
+  fRecTrackHist2(0)\r
+{\r
+  //\r
+  // constructor\r
+  //\r
+  for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { \r
+    fMCTrackHist1[i]=0;     \r
+    fMCPrimTrackHist1[i]=0;     \r
+    fMCPrimTrackHist2[i]=0;     \r
+    fMCSecTrackHist1[i]=0;     \r
+    fRecTrackHist1[i]=0;     \r
+    fRecTrackMultHist1[i]=0; \r
+  }\r
+\r
+  Init();\r
+}\r
+\r
+//_____________________________________________________________________________\r
+AlidNdPtAnalysisPbPb::~AlidNdPtAnalysisPbPb() {\r
+  //\r
+  // destructor\r
+  //\r
+  if(fTrackPtCorrelationMatrix) delete fTrackPtCorrelationMatrix; fTrackPtCorrelationMatrix=0;\r
+  //\r
+  if(fGenEventMatrix) delete fGenEventMatrix; fGenEventMatrix=0;\r
+  if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;\r
+  if(fRecEventMatrix) delete fRecEventMatrix; fRecEventMatrix=0;\r
+\r
+  //\r
+  if(fGenTrackEventMatrix) delete fGenTrackEventMatrix; fGenTrackEventMatrix=0;\r
+  if(fTriggerEventMatrix) delete fTriggerEventMatrix; fTriggerEventMatrix=0;\r
+  if(fRecTrackEventMatrix) delete fRecTrackEventMatrix; fRecTrackEventMatrix=0;\r
+\r
+  //\r
+  if(fGenTrackMatrix) delete fGenTrackMatrix; fGenTrackMatrix=0;\r
+  if(fGenPrimTrackMatrix) delete fGenPrimTrackMatrix; fGenPrimTrackMatrix=0;\r
+  if(fRecPrimTrackMatrix) delete fRecPrimTrackMatrix; fRecPrimTrackMatrix=0;\r
+  //\r
+  if(fRecTrackMatrix) delete fRecTrackMatrix; fRecTrackMatrix=0;\r
+  if(fRecSecTrackMatrix) delete fRecSecTrackMatrix; fRecSecTrackMatrix=0;\r
+  // \r
+  if(fRecMultTrackMatrix) delete fRecMultTrackMatrix; fRecMultTrackMatrix=0;\r
+  //\r
+  // Control histograms\r
+  //\r
+  if(fMCEventHist1) delete fMCEventHist1; fMCEventHist1=0;\r
+  if(fRecEventHist1) delete fRecEventHist1; fRecEventHist1=0;\r
+  if(fRecEventHist2) delete fRecEventHist2; fRecEventHist2=0;\r
+  if(fRecMCEventHist1) delete fRecMCEventHist1; fRecMCEventHist1=0;\r
+  if(fRecMCEventHist2) delete fRecMCEventHist2; fRecMCEventHist2=0;\r
+\r
+  //\r
+  for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { \r
+    if(fMCTrackHist1[i]) delete fMCTrackHist1[i]; fMCTrackHist1[i]=0;\r
+    if(fMCPrimTrackHist1[i]) delete fMCPrimTrackHist1[i]; fMCPrimTrackHist1[i]=0;\r
+    if(fMCPrimTrackHist2[i]) delete fMCPrimTrackHist2[i]; fMCPrimTrackHist2[i]=0;\r
+    if(fMCSecTrackHist1[i]) delete fMCSecTrackHist1[i]; fMCSecTrackHist1[i]=0;\r
+    if(fRecTrackHist1[i]) delete fRecTrackHist1[i]; fRecTrackHist1[i]=0;\r
+    if(fRecTrackMultHist1[i]) delete fRecTrackMultHist1[i]; fRecTrackMultHist1[i]=0;\r
+  }\r
+  if(fRecMCTrackHist1) delete fRecMCTrackHist1; fRecMCTrackHist1=0;\r
+  if(fMCMultRecTrackHist1) delete fMCMultRecTrackHist1; fMCMultRecTrackHist1=0; \r
+  if(fRecTrackHist2) delete fRecTrackHist2; fRecTrackHist2=0; \r
+  //\r
+  if(fAnalysisFolder) delete fAnalysisFolder; fAnalysisFolder=0;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AlidNdPtAnalysisPbPb::Init(){\r
+  //\r
+  // Init histograms\r
+  //\r
+\r
+  const Int_t multNbins = 1000;\r
+  const Int_t ptNbinsTrackEventCorr = 38;\r
+  const Int_t ptNbins = 55;\r
+  const Int_t etaNbins = 30;\r
+  const Int_t zvNbins = 12;\r
+\r
+  /*Double_t binsMult[multNbins+1] = {-0.5, 0.5 , 1.5 , 2.5 , 3.5 , 4.5 , 5.5 , 6.5 , 7.5 , 8.5,\r
+                                     9.5, 10.5, 11.5, 12.5, 13.5, 14.5, 15.5, 16.5, 17.5, 18.5,\r
+                                    19.5,20.5, 21.5, 22.5, 23.5, 24.5, 29.5, 149.5};\r
+  */\r
+\r
+  Double_t binsPtTrackEventCorr[ptNbinsTrackEventCorr+1] = {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,6.0,10.0,16.0};\r
+\r
+  Double_t binsPt[ptNbins+1] = {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};\r
+\r
+  Double_t binsEta[etaNbins+1] = {-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};\r
+\r
+  Double_t binsZv[zvNbins+1] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};\r
+  //Double_t binsMult[multNbins+1] = {0.,1.,2.,3.,4.,5.,6.,7.,8.,9.,10.,12.,14.,16.,18.,20.,30.,40.,50.,70.,90.,110.,150.};\r
+\r
+  //Int_t binsTrackMatrix[3]={zvNbins,ptNbins,etaNbins};\r
+  Int_t binsTrackEventCorrMatrix[3]={zvNbins,ptNbinsTrackEventCorr,etaNbins};\r
+\r
+  Int_t binsTrackPtCorrelationMatrix[3]={ptNbins,ptNbins,etaNbins};\r
+  fTrackPtCorrelationMatrix = new THnSparseF("fTrackPtCorrelationMatrix","Pt:mcPt:mcEta",3,binsTrackPtCorrelationMatrix);\r
+  fTrackPtCorrelationMatrix->SetBinEdges(0,binsPt);\r
+  fTrackPtCorrelationMatrix->SetBinEdges(1,binsPt);\r
+  fTrackPtCorrelationMatrix->SetBinEdges(2,binsEta);\r
+  fTrackPtCorrelationMatrix->GetAxis(0)->SetTitle("Pt (GeV/c)");\r
+  fTrackPtCorrelationMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fTrackPtCorrelationMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fTrackPtCorrelationMatrix->Sumw2();\r
+\r
+  //\r
+  // Efficiency and contamination correction matrices\r
+  //\r
+  Int_t binsEventMatrix[2]={zvNbins,multNbins};\r
+  Double_t minEventMatrix[2]={-25.,-0.5}; \r
+  Double_t maxEventMatrix[2]={25.,multNbins-0.5 }; \r
+\r
+  fGenEventMatrix = new THnSparseF("fGenEventMatrix","mcZv:multMB",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
+  fGenEventMatrix->SetBinEdges(0,binsZv);\r
+//  fGenEventMatrix->SetBinEdges(1,binsMult);\r
+  fGenEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fGenEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");\r
+  fGenEventMatrix->Sumw2();\r
+  //\r
+  fTriggerEventMatrix = new THnSparseF("fTriggerEventMatrix","mcZv:multMB",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
+  fTriggerEventMatrix->SetBinEdges(0,binsZv);\r
+//  fTriggerEventMatrix->SetBinEdges(1,binsMult);\r
+  fTriggerEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fTriggerEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");\r
+  fTriggerEventMatrix->Sumw2();\r
+  //\r
+  fRecEventMatrix = new THnSparseF("fRecEventMatrix","mcZv:multMB",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
+  fRecEventMatrix->SetBinEdges(0,binsZv);\r
+//  fRecEventMatrix->SetBinEdges(1,binsMult);\r
+  fRecEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fRecEventMatrix->GetAxis(1)->SetTitle("multiplicity MB");\r
+  fRecEventMatrix->Sumw2();\r
+\r
+  // \r
+  // track to event corrections\r
+  //\r
+  fGenTrackEventMatrix = new THnSparseF("fGenTrackEventMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fGenTrackEventMatrix->SetBinEdges(0,binsZv);\r
+  fGenTrackEventMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fGenTrackEventMatrix->SetBinEdges(2,binsEta);\r
+  fGenTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fGenTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fGenTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fGenTrackEventMatrix->Sumw2();\r
+  //\r
+  fTriggerTrackEventMatrix = new THnSparseF("fTriggerTrackEventMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fTriggerTrackEventMatrix->SetBinEdges(0,binsZv);\r
+  fTriggerTrackEventMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fTriggerTrackEventMatrix->SetBinEdges(2,binsEta);\r
+  fTriggerTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fTriggerTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fTriggerTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fTriggerTrackEventMatrix->Sumw2();\r
+  //\r
+  fRecTrackEventMatrix = new THnSparseF("fRecTrackEventMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fRecTrackEventMatrix->SetBinEdges(0,binsZv);\r
+  fRecTrackEventMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fRecTrackEventMatrix->SetBinEdges(2,binsEta);\r
+  fRecTrackEventMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fRecTrackEventMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fRecTrackEventMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fRecTrackEventMatrix->Sumw2();\r
+  //\r
+  // tracks correction matrices\r
+  //\r
+\r
+  //fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);\r
+  fGenTrackMatrix = new THnSparseF("fGenTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fGenTrackMatrix->SetBinEdges(0,binsZv);\r
+  //fGenTrackMatrix->SetBinEdges(1,binsPt);\r
+  fGenTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fGenTrackMatrix->SetBinEdges(2,binsEta);\r
+  fGenTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fGenTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fGenTrackMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fGenTrackMatrix->Sumw2();\r
+\r
+  //fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);\r
+  fGenPrimTrackMatrix = new THnSparseF("fGenPrimTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fGenPrimTrackMatrix->SetBinEdges(0,binsZv);\r
+  //fGenPrimTrackMatrix->SetBinEdges(1,binsPt);\r
+  fGenPrimTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fGenPrimTrackMatrix->SetBinEdges(2,binsEta);\r
+  fGenPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fGenPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fGenPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fGenPrimTrackMatrix->Sumw2();\r
+\r
+  //fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);\r
+  fRecPrimTrackMatrix = new THnSparseF("fRecPrimTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fRecPrimTrackMatrix->SetBinEdges(0,binsZv);\r
+  //fRecPrimTrackMatrix->SetBinEdges(1,binsPt);\r
+  fRecPrimTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fRecPrimTrackMatrix->SetBinEdges(2,binsEta);\r
+  fRecPrimTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fRecPrimTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fRecPrimTrackMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fRecPrimTrackMatrix->Sumw2();\r
+\r
+  //\r
+  //fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);\r
+  fRecTrackMatrix = new THnSparseF("fRecTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fRecTrackMatrix->SetBinEdges(0,binsZv);\r
+  //fRecTrackMatrix->SetBinEdges(1,binsPt);\r
+  fRecTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fRecTrackMatrix->SetBinEdges(2,binsEta);\r
+  fRecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fRecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)");\r
+  fRecTrackMatrix->GetAxis(2)->SetTitle("Eta");\r
+  fRecTrackMatrix->Sumw2();\r
+\r
+  //fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);\r
+  fRecSecTrackMatrix = new THnSparseF("fRecSecTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fRecSecTrackMatrix->SetBinEdges(0,binsZv);\r
+  //fRecSecTrackMatrix->SetBinEdges(1,binsPt);\r
+  fRecSecTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fRecSecTrackMatrix->SetBinEdges(2,binsEta);\r
+  fRecSecTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fRecSecTrackMatrix->GetAxis(1)->SetTitle("Pt (GeV/c)");\r
+  fRecSecTrackMatrix->GetAxis(2)->SetTitle("Eta");\r
+  fRecSecTrackMatrix->Sumw2();\r
+\r
+  //\r
+  //fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackMatrix);\r
+  fRecMultTrackMatrix = new THnSparseF("fRecMultTrackMatrix","mcZv:mcPt:mcEta",3,binsTrackEventCorrMatrix);\r
+  fRecMultTrackMatrix->SetBinEdges(0,binsZv);\r
+  //fRecMultTrackMatrix->SetBinEdges(1,binsPt);\r
+  fRecMultTrackMatrix->SetBinEdges(1,binsPtTrackEventCorr);\r
+  fRecMultTrackMatrix->SetBinEdges(2,binsEta);\r
+  fRecMultTrackMatrix->GetAxis(0)->SetTitle("mcZv (cm)");\r
+  fRecMultTrackMatrix->GetAxis(1)->SetTitle("mcPt (GeV/c)");\r
+  fRecMultTrackMatrix->GetAxis(2)->SetTitle("mcEta");\r
+  fRecMultTrackMatrix->Sumw2();\r
+\r
+  //\r
+  // Control analysis histograms\r
+  //\r
+\r
+  Int_t binsMCEventHist1[3]={100,100,140};\r
+  Double_t minMCEventHist1[3]={-0.1,-0.1,-35.}; \r
+  Double_t maxMCEventHist1[3]={0.1,0.1,35.}; \r
+  fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv",3,binsMCEventHist1,minMCEventHist1,maxMCEventHist1);\r
+  fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)");\r
+  fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)");\r
+  fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)");\r
+  fMCEventHist1->Sumw2();\r
+\r
+  //\r
+  Int_t binsRecEventHist1[3]={100,100,140};\r
+  Double_t minRecEventHist1[3]={-3.,-3.,-35.}; \r
+  Double_t maxRecEventHist1[3]={3.,3.,35.}; \r
+  \r
+  fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv",3,binsRecEventHist1,minRecEventHist1,maxRecEventHist1);\r
+  fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)");\r
+  fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)");\r
+  fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)");\r
+  fRecEventHist1->Sumw2();\r
+\r
+  //\r
+  Int_t binsRecEventHist2[3]={zvNbins,multNbins,multNbins};\r
+  Double_t minRecEventHist2[3]={-25.,-0.5,-0.5}; \r
+  Double_t maxRecEventHist2[3]={25.,multNbins-0.5,multNbins-0.5}; \r
+  \r
+  fRecEventHist2 = new THnSparseF("fRecEventHist2","Zv:multMB:mult",3,binsRecEventHist2,minRecEventHist2,maxRecEventHist2);\r
+  fRecEventHist2->SetBinEdges(0,binsZv);\r
+  fRecEventHist2->GetAxis(0)->SetTitle("Zv (cm)");\r
+  fRecEventHist2->GetAxis(1)->SetTitle("multiplicity MB");\r
+  fRecEventHist2->GetAxis(2)->SetTitle("multiplicity");\r
+  fRecEventHist2->Sumw2();\r
+\r
+  //\r
+  Double_t kFact = 0.1;\r
+  Int_t binsRecMCEventHist1[3]={100,100,100};\r
+  Double_t minRecMCEventHist1[3]={-10.0*kFact,-10.0*kFact,-10.0*kFact}; \r
+  Double_t maxRecMCEventHist1[3]={10.0*kFact,10.0*kFact,10.0*kFact}; \r
+   \r
+  fRecMCEventHist1 = new THnSparseF("fRecMCEventHist1","Xv-mcXv:Yv-mcYv:Zv-mcZv",3,binsRecMCEventHist1,minRecMCEventHist1,maxRecMCEventHist1);\r
+  fRecMCEventHist1->GetAxis(0)->SetTitle("Xv-mcXv (cm)");\r
+  fRecMCEventHist1->GetAxis(1)->SetTitle("Yv-mcYv (cm)");\r
+  fRecMCEventHist1->GetAxis(2)->SetTitle("Zv-mcZv (cm)");\r
+  fRecMCEventHist1->Sumw2();\r
+\r
+  //\r
+  Int_t binsRecMCEventHist2[3]={100,100,multNbins};\r
+  Double_t minRecMCEventHist2[3]={-10.0*kFact,-10.0*kFact,-0.5}; \r
+  Double_t maxRecMCEventHist2[3]={10.0*kFact,10.0*kFact,multNbins-0.5}; \r
+\r
+  fRecMCEventHist2 = new THnSparseF("fRecMCEventHist2","Xv-mcXv:Zv-mcZv:mult",3,binsRecMCEventHist2,minRecMCEventHist2,maxRecMCEventHist2);\r
+  fRecMCEventHist2->GetAxis(0)->SetTitle("Xv-mcXv (cm)");\r
+  fRecMCEventHist2->GetAxis(1)->SetTitle("Zv-mcZv (cm)");\r
+  fRecMCEventHist2->GetAxis(2)->SetTitle("multiplicity");\r
+  fRecMCEventHist2->Sumw2();\r
+  //\r
+  char name[256];\r
+  char title[256];\r
+  for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) \r
+  {\r
+  // THnSparse track histograms\r
\r
+  Int_t binsMCTrackHist1[3]=  {ptNbins, etaNbins, 90};\r
+  Double_t minMCTrackHist1[3]={0.,-1.,0.}; \r
+  Double_t maxMCTrackHist1[3]={10.,1.,2.*TMath::Pi()}; \r
+  sprintf(name,"fMCTrackHist1_%d",i);\r
+  sprintf(title,"mcPt:mcEta:mcPhi");\r
+  \r
+  fMCTrackHist1[i] = new THnSparseF(name,title,3,binsMCTrackHist1,minMCTrackHist1,maxMCTrackHist1);\r
+  fMCTrackHist1[i]->SetBinEdges(0,binsPt);\r
+  fMCTrackHist1[i]->SetBinEdges(1,binsEta);\r
+  fMCTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
+  fMCTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");\r
+  fMCTrackHist1[i]->GetAxis(2)->SetTitle("mcPhi (rad)");\r
+  fMCTrackHist1[i]->Sumw2();\r
+\r
+  Int_t binsMCPrimTrackHist1[5]=  {ptNbins,etaNbins,6,20,4000};\r
+  Double_t minMCPrimTrackHist1[5]={0.,-1.,0.,0.,0.}; \r
+  Double_t maxMCPrimTrackHist1[5]={10.,1.,6.,20.,4000.}; \r
+  sprintf(name,"fMCPrimTrackHist1_%d",i);\r
+  sprintf(title,"mcPt:mcEta:pid:mech:mother");\r
+  \r
+  fMCPrimTrackHist1[i] = new THnSparseF(name,title,5,binsMCPrimTrackHist1,minMCPrimTrackHist1,maxMCPrimTrackHist1);\r
+  fMCPrimTrackHist1[i]->SetBinEdges(0,binsPt);\r
+  fMCPrimTrackHist1[i]->SetBinEdges(1,binsEta);\r
+  fMCPrimTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
+  fMCPrimTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");\r
+  fMCPrimTrackHist1[i]->GetAxis(2)->SetTitle("pid");\r
+  fMCPrimTrackHist1[i]->GetAxis(3)->SetTitle("mech");\r
+  fMCPrimTrackHist1[i]->GetAxis(4)->SetTitle("mother");\r
+  fMCPrimTrackHist1[i]->Sumw2();\r
+\r
+  Int_t binsMCPrimTrackHist2[5]=  {4000,20,4000};\r
+  Double_t minMCPrimTrackHist2[5]={0.,0.,0.}; \r
+  Double_t maxMCPrimTrackHist2[5]={4000.,20.,4000.}; \r
+  sprintf(name,"fMCPrimTrackHist2_%d",i);\r
+  sprintf(title,"pdg:mech:mother");\r
+  \r
+  fMCPrimTrackHist2[i] = new THnSparseF(name,title,5,binsMCPrimTrackHist2,minMCPrimTrackHist2,maxMCPrimTrackHist2);\r
+  fMCPrimTrackHist2[i]->GetAxis(0)->SetTitle("pdg");\r
+  fMCPrimTrackHist2[i]->GetAxis(1)->SetTitle("mech");\r
+  fMCPrimTrackHist2[i]->GetAxis(2)->SetTitle("mother");\r
+  fMCPrimTrackHist2[i]->Sumw2();\r
+\r
+  Int_t binsMCSecTrackHist1[5]=  {ptNbins,etaNbins,6,20,4000};\r
+  Double_t minMCSecTrackHist1[5]={0.,-1.,0.,0.,0.}; \r
+  Double_t maxMCSecTrackHist1[5]={10.,1.,6.,20.,4000.}; \r
+  sprintf(name,"fMCSecTrackHist1_%d",i);\r
+  sprintf(title,"mcPt:mcEta:mcPhi:pid:mech:mother");\r
+  \r
+  fMCSecTrackHist1[i] = new THnSparseF(name,title,5,binsMCSecTrackHist1,minMCSecTrackHist1,maxMCSecTrackHist1);\r
+  fMCSecTrackHist1[i]->SetBinEdges(0,binsPt);\r
+  fMCSecTrackHist1[i]->SetBinEdges(1,binsEta);\r
+  fMCSecTrackHist1[i]->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
+  fMCSecTrackHist1[i]->GetAxis(1)->SetTitle("mcEta");\r
+  fMCSecTrackHist1[i]->GetAxis(2)->SetTitle("pid");\r
+  fMCSecTrackHist1[i]->GetAxis(3)->SetTitle("mech");\r
+  fMCSecTrackHist1[i]->GetAxis(4)->SetTitle("mother");\r
+  fMCSecTrackHist1[i]->Sumw2();\r
+\r
+  //\r
+\r
+  // \r
+  \r
+  Int_t binsRecTrackHist1[3]={ptNbins,etaNbins,90};\r
+  Double_t minRecTrackHist1[3]={0.,-1.,0.}; \r
+  Double_t maxRecTrackHist1[3]={10.,1.,2.*TMath::Pi()};\r
+  sprintf(name,"fRecTrackHist1_%d",i);\r
+  sprintf(title,"Pt:Eta:Phi");\r
+  fRecTrackHist1[i] = new THnSparseF(name,title,3,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1);\r
+  fRecTrackHist1[i]->SetBinEdges(0,binsPt);\r
+  fRecTrackHist1[i]->SetBinEdges(1,binsEta);\r
+  fRecTrackHist1[i]->GetAxis(0)->SetTitle("p_{T} (GeV/c)");\r
+  fRecTrackHist1[i]->GetAxis(1)->SetTitle("#eta");\r
+  fRecTrackHist1[i]->GetAxis(2)->SetTitle("#phi (rad)");\r
+  fRecTrackHist1[i]->Sumw2();\r
+\r
+  // \r
+  Int_t binsRecTrackMultHist1[2]={ptNbins,multNbins};\r
+  Double_t minRecTrackMultHist1[2]={0.,-0.5}; \r
+  Double_t maxRecTrackMultHist1[2]={10.,multNbins-0.5};\r
+  sprintf(name,"fRecTrackMultHist_%d",i);\r
+  sprintf(title,"Pt:Mult");\r
+  fRecTrackMultHist1[i] = new THnSparseF(name,title,2,binsRecTrackMultHist1,minRecTrackMultHist1,maxRecTrackMultHist1);\r
+  fRecTrackMultHist1[i]->SetBinEdges(0,binsPt);\r
+  fRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");\r
+  fRecTrackMultHist1[i]->GetAxis(1)->SetTitle("multiplicity");\r
+  fRecTrackMultHist1[i]->Sumw2();\r
+  }\r
+\r
+  Int_t binsRecMCTrackHist1[4] = {ptNbins,etaNbins,100,100};\r
+  Double_t minRecMCTrackHist1[4]={0.,-1.,-0.5,-0.5}; \r
+  Double_t maxRecMCTrackHist1[4]={20.,1.,0.5,0.5}; \r
+  sprintf(name,"fRecMCTrackHist1");\r
+  sprintf(title,"mcPt:mcEta:(Pt-mcPt)/mcPt:(Eta-mcEta)");\r
+  fRecMCTrackHist1 = new THnSparseF(name,title,4,binsRecMCTrackHist1,minRecMCTrackHist1,maxRecMCTrackHist1);\r
+  fRecMCTrackHist1->SetBinEdges(0,binsPt);\r
+  fRecMCTrackHist1->SetBinEdges(1,binsEta);\r
+  fRecMCTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
+  fRecMCTrackHist1->GetAxis(1)->SetTitle("mcEta");\r
+  fRecMCTrackHist1->GetAxis(2)->SetTitle("(Pt-mcPt)/mcPt");\r
+  fRecMCTrackHist1->GetAxis(3)->SetTitle("Eta-mcEta");\r
+\r
+  Int_t binsMCMultRecTrackHist1[3] = {ptNbins,etaNbins,6};\r
+  Double_t minMCMultRecTrackHist1[3]={0.,-1.,0.}; \r
+  Double_t maxMCMultRecTrackHist1[3]={20.,1.,6.}; \r
+  sprintf(name,"fMCMultRecTrackHist1");\r
+  sprintf(title,"mcPt:mcEta:pid");\r
+  fMCMultRecTrackHist1 = new THnSparseF(name,title,3,binsMCMultRecTrackHist1,minMCMultRecTrackHist1,maxMCMultRecTrackHist1);\r
+  fMCMultRecTrackHist1->SetBinEdges(0,binsPt);\r
+  fMCMultRecTrackHist1->SetBinEdges(1,binsEta);\r
+  fMCMultRecTrackHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
+  fMCMultRecTrackHist1->GetAxis(1)->SetTitle("mcEta");\r
+  fMCMultRecTrackHist1->GetAxis(2)->SetTitle("pid");\r
+\r
+  //nClust:chi2PerClust:pt:eta:phi\r
+  Int_t binsRecTrackHist2[5]={160,100,ptNbins,etaNbins,90};\r
+  Double_t minRecTrackHist2[5]={0., 0., 0, -1.5, 0.};\r
+  Double_t maxRecRecTrackHist2[5]={160.,10., 16, 1.5, 2.*TMath::Pi()};\r
+\r
+  fRecTrackHist2 = new THnSparseF("fRecTrackHist2","nClust:chi2PerClust:pt:eta:phi",5,binsRecTrackHist2,minRecTrackHist2,maxRecRecTrackHist2);\r
+  fRecTrackHist2->SetBinEdges(2,binsPt);\r
+  fRecTrackHist2->SetBinEdges(3,binsEta);\r
+  fRecTrackHist2->GetAxis(0)->SetTitle("nClust");\r
+  fRecTrackHist2->GetAxis(1)->SetTitle("chi2PerClust");\r
+  fRecTrackHist2->GetAxis(2)->SetTitle("p_{T} (GeV/c)");\r
+  fRecTrackHist2->GetAxis(3)->SetTitle("#eta");\r
+  fRecTrackHist2->GetAxis(4)->SetTitle("#phi (rad)");\r
+  fRecTrackHist2->Sumw2();\r
+\r
+  // init folder\r
+  fAnalysisFolder = CreateFolder("folderdNdPt","Analysis dNdPt Folder");\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AlidNdPtAnalysisPbPb::Process(AliESDEvent *const esdEvent, AliMCEvent *const mcEvent)\r
+{\r
+  //\r
+  // Process real and/or simulated events\r
+  //\r
+  if(!esdEvent) {\r
+    AliDebug(AliLog::kError, "esdEvent not available");\r
+    return;\r
+  }\r
+\r
+\r
+  // get selection cuts\r
+  AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
+  AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
+  AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
+\r
+  if(!evtCuts || !accCuts  || !esdTrackCuts) {\r
+    AliDebug(AliLog::kError, "cuts not available");\r
+    return;\r
+  }\r
+\r
+  // trigger selection\r
+  Bool_t isEventTriggered = kTRUE;\r
+  AliPhysicsSelection *trigSel = NULL;\r
+  AliTriggerAnalysis *trigAna = NULL;\r
+\r
+  if(evtCuts->IsTriggerRequired())  \r
+  {\r
+    //\r
+    trigSel = GetPhysicsTriggerSelection();\r
+    if(!trigSel) {\r
+      printf("cannot get trigSel \n");\r
+      return;\r
+    }\r
+\r
+    //\r
+    if(IsUseMCInfo()) \r
+    { \r
+      trigSel->SetAnalyzeMC();\r
+\r
+      \r
+        isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);\r
+       \r
+        if(GetTrigger() == AliTriggerAnalysis::kV0AND) \r
+       {\r
+          trigAna = trigSel->GetTriggerAnalysis();\r
+          if(!trigAna) \r
+            return;\r
+\r
+          isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());\r
+        }//if(GetTrigger() == AliTriggerAnalysis::kV0AND)\r
+     }//if(IsUseMCInfo())\r
+  }//if(evtCuts->IsTriggerRequired())\r
+\r
+\r
+\r
+  // use MC information\r
+  AliHeader* header = 0;\r
+  AliGenEventHeader* genHeader = 0;\r
+  AliStack* stack = 0;\r
+  TArrayF vtxMC(3);\r
+\r
+  Int_t multMCTrueTracks = 0;\r
+  if(IsUseMCInfo())\r
+  {\r
+    //\r
+    if(!mcEvent) {\r
+      AliDebug(AliLog::kError, "mcEvent not available");\r
+      return;\r
+    }\r
+    // get MC event header\r
+    header = mcEvent->Header();\r
+    if (!header) {\r
+      AliDebug(AliLog::kError, "Header not available");\r
+      return;\r
+    }\r
+    // MC particle stack\r
+    stack = mcEvent->Stack();\r
+    if (!stack) {\r
+      AliDebug(AliLog::kError, "Stack not available");\r
+      return;\r
+    }\r
+\r
+    // get MC vertex\r
+    genHeader = header->GenEventHeader();\r
+    if (!genHeader) {\r
+      AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");\r
+      return;\r
+    }\r
+    genHeader->PrimaryVertex(vtxMC);\r
+\r
+    Double_t vMCEventHist1[3]={vtxMC[0],vtxMC[1],vtxMC[2]};\r
+    fMCEventHist1->Fill(vMCEventHist1);\r
+\r
+    // multipliticy of all MC primary tracks\r
+    // in Zv, pt and eta ranges)\r
+    multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);\r
+\r
+  } // end bUseMC\r
+\r
+  // get reconstructed vertex  \r
+  const AliESDVertex* vtxESD = 0; \r
+  Bool_t isRecVertex = kTRUE;\r
+  if(evtCuts->IsRecVertexRequired()) \r
+  {\r
+//    Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();\r
+    //Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();\r
+    //vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints); \r
+    //isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, esdEvent->GetPrimaryVertexSPD(), GetAnalysisMode(), kFALSE);\r
+    if(GetAnalysisMode() == AlidNdPtHelper::kTPC) {\r
+      vtxESD = esdEvent->GetPrimaryVertexTPC();\r
+      //isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, esdEvent->GetPrimaryVertexSPD(), GetAnalysisMode(), kFALSE);\r
+    }\r
+  }\r
+\r
+  Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex; \r
+  //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);\r
+  //printf("GetAnalysisMode() %d \n",GetAnalysisMode());\r
+\r
+  // vertex contributors\r
+  Int_t multMBTracks = 0; \r
+  if(GetAnalysisMode() == AlidNdPtHelper::kTPC) \r
+  {\r
+     if(vtxESD->GetStatus() && isRecVertex) {\r
+         //multMBTracks = AlidNdPtHelper::GetTPCMBTrackMult(esdEvent,evtCuts,accCuts,esdTrackCuts);\r
+         multMBTracks = vtxESD->GetNContributors();\r
+     }\r
+  } \r
+  else {\r
+    AliDebug(AliLog::kError, Form("Found analysis type %d", GetAnalysisMode()));\r
+    return; \r
+  }\r
+\r
+  TObjArray *allChargedTracks=0;\r
+  //Int_t multAll=0, multAcc=0, multRec=0;\r
+  Int_t multAll=0, multRec=0;\r
+  Int_t *labelsAll=0, *labelsAcc=0, *labelsRec=0;\r
+\r
+  // check event cuts\r
+  if(isEventOK && isEventTriggered)\r
+  {\r
+\r
+    // get all charged tracks\r
+    allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());\r
+    if(!allChargedTracks) return;\r
+\r
+    Int_t entries = allChargedTracks->GetEntries();\r
+\r
+    labelsAll = new Int_t[entries];\r
+    labelsAcc = new Int_t[entries];\r
+    labelsRec = new Int_t[entries];\r
+    for(Int_t i=0; i<entries;++i) \r
+    {\r
+      AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
+\r
+\r
+      if(!track) continue;\r
+      if(track->Charge()==0) continue;\r
+\r
+\r
+\r
+      // only postive charged \r
+      if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0) \r
+        continue;\r
+      \r
+      // only negative charged \r
+      if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0) \r
+        continue;\r
+\r
+      //\r
+      FillHistograms(track,stack,AlidNdPtHelper::kAllTracks); \r
+      labelsAll[multAll] = TMath::Abs(track->GetLabel());\r
+\r
+      multAll++;\r
+\r
+      if(esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track)) {\r
+\r
+         FillHistograms(track,stack,AlidNdPtHelper::kRecTracks); \r
+         labelsRec[multRec] = TMath::Abs(track->GetLabel());\r
+\r
+         multRec++;\r
+\r
+      }\r
+     }//loop over entries\r
+\r
+\r
+     // fill track multiplicity histograms\r
+     //FillHistograms(allChargedTracks,labelsAll,multAll,labelsAcc,multAcc,labelsRec,multRec);\r
+\r
+     Double_t vRecEventHist1[3] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv()};\r
+     fRecEventHist1->Fill(vRecEventHist1);\r
+\r
+     Double_t vRecEventHist2[3] = {vtxESD->GetZv(),multMBTracks,multRec};\r
+     fRecEventHist2->Fill(vRecEventHist2);\r
+\r
+   } // triggered and event vertex\r
+\r
+   if(IsUseMCInfo())  \r
+   {\r
+\r
+     // \r
+     // event level corrections (zv,N_MB)\r
+     //\r
+     // all inelastic\r
+     Double_t vEventMatrix[2] = {vtxMC[2],multMBTracks};\r
+     fGenEventMatrix->Fill(vEventMatrix); \r
+     if(isEventTriggered) fTriggerEventMatrix->Fill(vEventMatrix);\r
+     if(isEventOK && isEventTriggered) fRecEventMatrix->Fill(vEventMatrix);\r
+\r
+     //\r
+     // track-event level corrections (zv,pt,eta)\r
+     //\r
+     for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc) \r
+     {\r
+       TParticle* particle = stack->Particle(iMc);\r
+       if (!particle)\r
+       continue;\r
+\r
+       // only charged particles\r
+       if(!particle->GetPDG()) continue;\r
+       Double_t charge = particle->GetPDG()->Charge()/3.;\r
+       if ( TMath::Abs(charge) < 0.001 )\r
+        continue;\r
+\r
+       // only postive charged \r
+       if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.) \r
+        continue;\r
+       \r
+       // only negative charged \r
+       if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.) \r
+       continue;\r
+      \r
+       // physical primary\r
+       Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
+       if(!prim) continue;\r
+\r
+       // checked accepted\r
+       if(accCuts->AcceptTrack(particle)) \r
+       {\r
+         Double_t vTrackEventMatrix[3] = {vtxMC[2], particle->Pt(), particle->Eta()}; \r
+         fGenTrackEventMatrix->Fill(vTrackEventMatrix);\r
+\r
+       }\r
+     }//loop over stack\r
+\r
+     // \r
+     // track-level corrections (zv,pt,eta)\r
+     //\r
+     if(isEventOK && isEventTriggered)\r
+     {\r
+\r
+       // fill MC and rec event control histograms\r
+       if(fHistogramsOn) {\r
+         Double_t vRecMCEventHist1[3] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetYv()-vtxMC[1],vtxESD->GetZv()-vtxMC[2]};\r
+         fRecMCEventHist1->Fill(vRecMCEventHist1);//\r
+\r
+         Double_t vRecMCEventHist2[3] = {vtxESD->GetXv()-vtxMC[0],vtxESD->GetZv()-vtxMC[2],multMBTracks};\r
+         fRecMCEventHist2->Fill(vRecMCEventHist2);\r
+\r
+       }//\r
+\r
+       //\r
+       // MC histograms for track efficiency studies\r
+       //\r
+       for (Int_t iMc = 0; iMc < stack->GetNtrack(); ++iMc) \r
+       {\r
+         TParticle* particle = stack->Particle(iMc);\r
+         if (!particle)\r
+         continue;\r
+\r
+         Double_t vTrackMatrix[3] = {vtxMC[2],particle->Pt(),particle->Eta()}; \r
+\r
+         // only charged particles\r
+         if(!particle->GetPDG()) continue;\r
+         Double_t charge = particle->GetPDG()->Charge()/3.;\r
+         if (TMath::Abs(charge) < 0.001)\r
+         continue;\r
+\r
+         // only postive charged \r
+         if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.) \r
+        continue;\r
+       \r
+         // only negative charged \r
+         if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.) \r
+        continue;\r
+      \r
+         // physical primary\r
+         Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
+\r
+         // check accepted\r
+         if(accCuts->AcceptTrack(particle)) \r
+        {\r
+\r
+           if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) \r
+            fGenPrimTrackMatrix->Fill(vTrackMatrix);\r
+\r
+          // fill control histograms\r
+           if(fHistogramsOn) \r
+            FillHistograms(stack,iMc,AlidNdPtHelper::kAccTracks); \r
+\r
+           // check multiple found tracks\r
+          Int_t multCount = 0;\r
+           for(Int_t iRec=0; iRec<multRec; ++iRec)\r
+           {\r
+             if(iMc == labelsRec[iRec]) \r
+            {\r
+              multCount++;\r
+              if(multCount>1)\r
+              {  \r
+                 fRecMultTrackMatrix->Fill(vTrackMatrix);\r
+\r
+                 // fill control histogram\r
+                 if(fHistogramsOn) {\r
+                   Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);\r
+                   Double_t vMCMultRecTrackHist1[3] = {particle->Pt(), particle->Eta(), pid};\r
+                  fMCMultRecTrackHist1->Fill(vMCMultRecTrackHist1);\r
+                }\r
+              }\r
+            }\r
+          }\r
+\r
+           // check reconstructed\r
+           for(Int_t iRec=0; iRec<multRec; ++iRec)\r
+           {\r
+             if(iMc == labelsRec[iRec]) \r
+            {\r
+               fRecTrackMatrix->Fill(vTrackMatrix);\r
+\r
+               if( AlidNdPtHelper::IsPrimaryParticle(stack, iMc, GetParticleMode()) ) {\r
+                fRecPrimTrackMatrix->Fill(vTrackMatrix);\r
+              }\r
+               if(!prim) fRecSecTrackMatrix->Fill(vTrackMatrix);\r
+\r
+              // fill control histograms\r
+               if(fHistogramsOn) \r
+                 FillHistograms(stack,iMc,AlidNdPtHelper::kRecTracks); \r
+             \r
+               break;\r
+             }//if(iMc == labelsRec[iRec])\r
+           }//reco tracks\r
+         }//accepted tracks\r
+       }//stack loop\r
+     }//is triggered\r
+   } // end bUseMC\r
+\r
+  if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;\r
+  if(labelsAll) delete [] labelsAll; labelsAll = 0;\r
+  if(labelsAcc) delete [] labelsAcc; labelsAcc = 0;\r
+  if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
+\r
+  if(!evtCuts->IsRecVertexRequired() && vtxESD != NULL) delete vtxESD;\r
+\r
+}\r
+\r
+//_____________________________________________________________________________\r
+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) {\r
+ // multiplicity  histograms\r
\r
+\r
+  if(!allChargedTracks) return;\r
+  if(!labelsAll) return;\r
+  if(!labelsAcc) return;\r
+  if(!labelsRec) return;\r
+\r
+  Int_t entries = allChargedTracks->GetEntries();\r
+  for(Int_t i=0; i<entries; ++i)\r
+  {\r
+     AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
+     if(!track) continue;\r
+     if(track->Charge() == 0) continue;\r
+\r
+     Int_t label = TMath::Abs(track->GetLabel());\r
+     for(Int_t iAll=0; iAll<multAll; ++iAll) {\r
+       if(label == labelsAll[iAll]) {\r
+         Double_t v1[2] = {track->Pt(), multAll}; \r
+         fRecTrackMultHist1[AlidNdPtHelper::kAllTracks]->Fill(v1);\r
+       }\r
+     }\r
+     for(Int_t iAcc=0; iAcc<multAcc; ++iAcc) {\r
+       if(label == labelsAcc[iAcc]) {\r
+         Double_t v2[2] = {track->Pt(), multAcc}; \r
+         fRecTrackMultHist1[AlidNdPtHelper::kAccTracks]->Fill(v2);\r
+       }\r
+     }\r
+     for(Int_t iRec=0; iRec<multRec; ++iRec) {\r
+       if(label == labelsRec[iRec]) {\r
+         Double_t v3[2] = {track->Pt(), multRec}; \r
+         fRecTrackMultHist1[AlidNdPtHelper::kRecTracks]->Fill(v3);\r
+       }//out\r
+     }\r
+  }\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AlidNdPtAnalysisPbPb::FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, AlidNdPtHelper::TrackObject trackObj)\r
+{\r
+\r
+  //\r
+  // Fill ESD track and MC histograms \r
+  //\r
+  if(!esdTrack) return;\r
+\r
+  Float_t q = esdTrack->Charge();\r
+  if(TMath::Abs(q) < 0.001) return;\r
+\r
+  Float_t pt = esdTrack->Pt();\r
+  Float_t eta = esdTrack->Eta();\r
+  Float_t phi = esdTrack->Phi();\r
+\r
+  Float_t dca[2], bCov[3];\r
+  esdTrack->GetImpactParameters(dca,bCov);\r
+\r
+  Int_t nClust = esdTrack->GetTPCclusters(0);\r
+  Float_t chi2PerCluster = 0.;\r
+  if(nClust>0.) chi2PerCluster = esdTrack->GetTPCchi2()/Float_t(nClust);\r
+\r
+\r
+  // fill histograms\r
+  Double_t values[3] = {pt,eta,phi};     \r
+  fRecTrackHist1[trackObj]->Fill(values);\r
+\r
+  //\r
+  // Fill rec vs MC information\r
+  //\r
+  if(!stack) return;\r
+\r
+  Int_t label = TMath::Abs(esdTrack->GetLabel()); \r
+  //if(label == 0) return;\r
+\r
+  TParticle* particle = stack->Particle(label);\r
+  if(!particle) return;\r
+\r
+  Int_t motherPdg = -1;\r
+  TParticle* mother = 0;\r
+\r
+  //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);\r
+  Int_t motherLabel = particle->GetMother(0); \r
+  if(motherLabel>0) mother = stack->Particle(motherLabel);\r
+  if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only\r
+  //Int_t mech = particle->GetUniqueID(); // production mechanism\r
+\r
+  if(!particle->GetPDG()) return;\r
+  Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 \r
+  if(TMath::Abs(gq)<0.001) return;\r
+  Float_t gpt = particle->Pt();\r
+  Float_t geta = particle->Eta();\r
+\r
+  Double_t dpt=0;\r
+  //printf("pt %f, gpt %f \n",pt,gpt);\r
+  if(gpt) dpt = (pt-gpt)/gpt;\r
+  Double_t deta = (eta-geta);\r
\r
+  // fill histograms\r
+  if(trackObj == AlidNdPtHelper::kRecTracks)  //RecTracks???\r
+  {\r
+    Double_t vTrackPtCorrelationMatrix[3]={pt,gpt,geta};\r
+    fTrackPtCorrelationMatrix->Fill(vTrackPtCorrelationMatrix);\r
+\r
+    Double_t vRecMCTrackHist1[4]={gpt,geta,dpt,deta};\r
+    fRecMCTrackHist1->Fill(vRecMCTrackHist1);\r
+  }\r
+}\r
+\r
+//_____________________________________________________________________________\r
+void AlidNdPtAnalysisPbPb::FillHistograms(AliStack *const stack, Int_t label, AlidNdPtHelper::TrackObject trackObj)\r
+{\r
+\r
+  // Fill MC histograms\r
+  if(!stack) return;\r
+\r
+  TParticle* particle = stack->Particle(label);\r
+  if(!particle) return;\r
+\r
+  Int_t motherPdg = -1;\r
+  TParticle* mother = 0;\r
+\r
+  //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);\r
+  Int_t motherLabel = particle->GetMother(0); \r
+  if(motherLabel>0) mother = stack->Particle(motherLabel);\r
+  if(mother) motherPdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only\r
+  Int_t mech = particle->GetUniqueID(); // production mechanism\r
+\r
+  if(!particle->GetPDG()) return;\r
+  Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 \r
+  if(TMath::Abs(gq) < 0.001) return;\r
+\r
+  Float_t gpt = particle->Pt();\r
+  //Float_t qgpt = particle->Pt() * gq;\r
+  Float_t geta = particle->Eta();\r
+  Float_t gphi = particle->Phi();\r
+  //Float_t gpz = particle->Pz();\r
+\r
+  Bool_t prim = stack->IsPhysicalPrimary(label);\r
+  //Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();\r
+\r
+  Int_t pid = AlidNdPtHelper::ConvertPdgToPid(particle);\r
+\r
+  //if(prim&&pid==5) printf("pdgcode %d, production mech %d \n",particle->GetPdgCode(),mech);\r
+  //if(!prim) printf("motherPdg %d, particle %d, production mech %d\n",motherPdg, particle->GetPdgCode(),mech);\r
+  \r
+  //\r
+  // fill histogram\r
+  //\r
+  Double_t vMCTrackHist1[3] = {gpt,geta,gphi};\r
+  fMCTrackHist1[trackObj]->Fill(vMCTrackHist1);\r
+\r
+  Double_t vMCPrimTrackHist1[5] = {gpt,geta,pid,mech,motherPdg};\r
+  Double_t vMCPrimTrackHist2[5] = {TMath::Abs(particle->GetPdgCode()),mech,motherPdg};\r
+\r
+  //if(prim && AliPWG0Helper::IsPrimaryCharged(particle, label)) fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);\r
+\r
+  if(prim) { \r
+    fMCPrimTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);\r
+    if(pid == 5) fMCPrimTrackHist2[trackObj]->Fill(vMCPrimTrackHist2);\r
+  }\r
+  else { \r
+    fMCSecTrackHist1[trackObj]->Fill(vMCPrimTrackHist1);\r
+  }\r
+\r
+}\r
+\r
+//_____________________________________________________________________________\r
+Long64_t AlidNdPtAnalysisPbPb::Merge(TCollection* const list) \r
+{\r
+  // Merge list of objects (needed by PROOF)\r
+\r
+  if (!list)\r
+  return 0;\r
+\r
+  if (list->IsEmpty())\r
+  return 1;\r
+\r
+  TIterator* iter = list->MakeIterator();\r
+  TObject* obj = 0;\r
+\r
+  //\r
+  TList *collPhysSelection = new TList;\r
+\r
+  // collection of generated histograms\r
+\r
+  Int_t count=0;\r
+  while((obj = iter->Next()) != 0) {\r
+    AlidNdPtAnalysisPbPb* entry = dynamic_cast<AlidNdPtAnalysisPbPb*>(obj);\r
+    if (entry == 0) continue; \r
+\r
+    // physics selection\r
+    //printf("entry->GetPhysicsTriggerSelection() %p \n", entry->GetPhysicsTriggerSelection());\r
+    collPhysSelection->Add(entry->GetPhysicsTriggerSelection());\r
+    \r
+    //\r
+    fTrackPtCorrelationMatrix->Add(entry->fTrackPtCorrelationMatrix);\r
+\r
+    //\r
+    fGenEventMatrix->Add(entry->fGenEventMatrix);\r
+\r
+    fTriggerEventMatrix->Add(entry->fTriggerEventMatrix);\r
+\r
+    fRecEventMatrix->Add(entry->fRecEventMatrix);\r
+    //\r
+    fGenTrackEventMatrix->Add(entry->fGenTrackEventMatrix);\r
+\r
+    fTriggerTrackEventMatrix->Add(entry->fTriggerTrackEventMatrix);\r
+\r
+    fRecTrackEventMatrix->Add(entry->fRecTrackEventMatrix);\r
+\r
+    //\r
+    fGenTrackMatrix->Add(entry->fGenTrackMatrix);\r
+    fGenPrimTrackMatrix->Add(entry->fGenPrimTrackMatrix);\r
+    fRecPrimTrackMatrix->Add(entry->fRecPrimTrackMatrix);\r
+    //\r
+    fRecTrackMatrix->Add(entry->fRecTrackMatrix);\r
+    fRecSecTrackMatrix->Add(entry->fRecSecTrackMatrix);\r
+    //\r
+    fRecMultTrackMatrix->Add(entry->fRecMultTrackMatrix);\r
+\r
+    //\r
+    // control analysis histograms\r
+    //\r
+    fMCEventHist1->Add(entry->fMCEventHist1);\r
+    fRecEventHist1->Add(entry->fRecEventHist1);\r
+    fRecEventHist2->Add(entry->fRecEventHist2);\r
+    fRecMCEventHist1->Add(entry->fRecMCEventHist1);\r
+    fRecMCEventHist2->Add(entry->fRecMCEventHist2);\r
+\r
+\r
+    for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {\r
+      fMCTrackHist1[i]->Add(entry->fMCTrackHist1[i]);\r
+\r
+      fMCPrimTrackHist1[i]->Add(entry->fMCPrimTrackHist1[i]);\r
+      fMCPrimTrackHist2[i]->Add(entry->fMCPrimTrackHist2[i]);\r
+      fMCSecTrackHist1[i]->Add(entry->fMCSecTrackHist1[i]);\r
+\r
+      fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);\r
+      fRecTrackMultHist1[i]->Add(entry->fRecTrackMultHist1[i]);\r
+    }\r
+    fRecMCTrackHist1->Add(entry->fRecMCTrackHist1);\r
+    fMCMultRecTrackHist1->Add(entry->fMCMultRecTrackHist1);\r
+    fRecTrackHist2->Add(entry->fRecTrackHist2);\r
+\r
+  count++;\r
+  }\r
+\r
+  AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();\r
+  trigSelection->Merge(collPhysSelection);\r
+  if(collPhysSelection) delete collPhysSelection;\r
+\r
+return count;\r
+}\r
+\r
+\r
+\r
+//_____________________________________________________________________________\r
+void AlidNdPtAnalysisPbPb::Analyse() \r
+{\r
+  // Analyse histograms\r
+  //\r
+  TH1::AddDirectory(kFALSE);\r
+  TH1 *h=0, *h1=0, *h2=0, *h2c = 0; \r
+  THnSparse *hs=0; \r
+  TH2 *h2D=0; \r
+\r
+  char name[256];\r
+  TObjArray *aFolderObj = new TObjArray;\r
+  \r
+  //\r
+  // LHC backgraund in all and 0-bins\r
+  //\r
+  AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection();\r
+  trigSel->SaveHistograms("physics_selection");\r
+\r
+  //\r
+  // Reconstructed event vertex\r
+  //\r
+  h = fRecEventHist1->Projection(0);\r
+  h->SetName("Xv");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecEventHist1->Projection(1);\r
+  h->SetName("Yv");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecEventHist1->Projection(2);\r
+  h->SetName("Zv");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // multiplicity\r
+  //\r
+  h = fRecEventHist2->Projection(1);\r
+  h->SetName("multMB");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecEventHist2->Projection(2);\r
+  h->SetName("multiplicity");\r
+  aFolderObj->Add(h);\r
+\r
+  h2D = fRecEventHist2->Projection(0,1); \r
+  h2D->SetName("Zv_vs_multiplicity_MB");\r
+  aFolderObj->Add(h2D);\r
+\r
+  //\r
+  // reconstructed pt histograms\r
+  //\r
+  h = fRecTrackHist1[0]->Projection(0);\r
+  h->Scale(1.,"width");\r
+  h->SetName("pt_all_ch");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecTrackHist1[1]->Projection(0);\r
+  h->Scale(1.,"width");\r
+  h->SetName("pt_acc");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecTrackHist1[2]->Projection(0);\r
+  h->Scale(1.,"width");\r
+  h->SetName("pt_rec");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // reconstructed eta histograms\r
+  //\r
+  h = fRecTrackHist1[0]->Projection(1);\r
+  h->SetName("eta_all_ch");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecTrackHist1[1]->Projection(1);\r
+  h->SetName("eta_acc");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecTrackHist1[2]->Projection(1);\r
+  h->SetName("eta_rec");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // reconstructed phi histograms\r
+  //\r
+  h = fRecTrackHist1[0]->Projection(2);\r
+  h->SetName("phi_all_ch");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecTrackHist1[1]->Projection(2);\r
+  h->SetName("phi_acc");\r
+  aFolderObj->Add(h);\r
+\r
+  h = fRecTrackHist1[2]->Projection(2);\r
+  h->SetName("phi_rec");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // reconstructed eta:pt histograms\r
+  //\r
+  h2D = fRecTrackHist1[0]->Projection(1,0);\r
+  h2D->SetName("pt_eta_all_ch");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecTrackHist1[1]->Projection(1,0);\r
+  h2D->SetName("pt_eta_acc");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecTrackHist1[2]->Projection(1,0);\r
+  h2D->SetName("pt_eta_rec");\r
+  aFolderObj->Add(h2D);\r
+\r
+  //\r
+  // reconstructed phi:pt histograms\r
+  //\r
+  h2D = fRecTrackHist1[0]->Projection(2,0);\r
+  h2D->SetName("pt_phi_all_ch");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecTrackHist1[1]->Projection(2,0);\r
+  h2D->SetName("pt_phi_acc");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecTrackHist1[2]->Projection(2,0);\r
+  h2D->SetName("pt_phi_rec");\r
+  aFolderObj->Add(h2D);\r
+\r
+  //\r
+  // reconstructed phi:eta histograms\r
+  //\r
+  h2D = fRecTrackHist1[0]->Projection(2,1);\r
+  h2D->SetName("eta_phi_all_ch");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecTrackHist1[1]->Projection(2,1);\r
+  h2D->SetName("eta_phi_acc");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecTrackHist1[2]->Projection(2,1);\r
+  h2D->SetName("eta_phi_rec");\r
+  aFolderObj->Add(h2D);\r
+\r
+  //\r
+  // reconstructed nClust, chi2 vs pt, eta, phi\r
+  //\r
+  if(fHistogramsOn) {\r
+\r
+    h2D = fRecTrackHist2->Projection(0,1);\r
+    h2D->SetName("nClust_chi2_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+    h2D = fRecTrackHist2->Projection(0,2);\r
+    h2D->SetName("nClust_pt_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+    h2D = fRecTrackHist2->Projection(0,3);\r
+    h2D->SetName("nClust_eta_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+    h2D = fRecTrackHist2->Projection(0,4);\r
+    h2D->SetName("nClust_phi_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+    h2D = fRecTrackHist2->Projection(1,2);\r
+    h2D->SetName("chi2_pt_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+    h2D = fRecTrackHist2->Projection(1,3);\r
+    h2D->SetName("chi2_eta_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+    h2D = fRecTrackHist2->Projection(1,4);\r
+    h2D->SetName("chi2_phi_rec");\r
+    aFolderObj->Add(h2D);\r
+\r
+  }\r
+\r
+  //\r
+  // calculate corrections for empty events\r
+  // with multMB==0 \r
+  //\r
+\r
+  //\r
+  // normalised zv to generate zv for triggered events\r
+  //\r
+  h = fRecEventHist2->Projection(0);\r
+  if( h->Integral() ) h->Scale(1./h->Integral());\r
+  h->SetName("zv_distribution_norm");\r
+  aFolderObj->Add(h);\r
\r
+  //\r
+  // MC available\r
+  //\r
+  if(IsUseMCInfo()) {\r
+\r
+  //\r
+  // Event vertex resolution\r
+  //\r
+  h2D = fRecMCEventHist2->Projection(0,2);\r
+  h2D->SetName("DeltaXv_vs_mult");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = fRecMCEventHist2->Projection(1,2);\r
+  h2D->SetName("DeltaZv_vs_mult");\r
+  aFolderObj->Add(h2D);\r
+\r
+  //\r
+  // normalised zv to get trigger/trigger+vertex event differences\r
+  // F(zv) = E_trig(zv,0)/Int(E_trig(zv,0) / Sum(E_trigvtx(zv,n))/Sum(Int(E_trigvtx(zv,n))dzv)\r
+  //\r
+  fTriggerEventMatrix->GetAxis(1)->SetRangeUser(0.,0.);\r
+  h = fTriggerEventMatrix->Projection(0);\r
+  h2D = fTriggerEventMatrix->Projection(0,1);\r
+  if(h2D->Integral()) h->Scale(1./h2D->Integral());\r
+\r
+  h1 = fRecEventMatrix->Projection(0);\r
+  h2D = fRecEventMatrix->Projection(0,1);\r
+  if(h2D->Integral()) h1->Scale(1./h2D->Integral());\r
+\r
+  h->Divide(h1);\r
+  h->SetName("zv_empty_events_norm");\r
+  aFolderObj->Add(h);\r
+  \r
+  fTriggerEventMatrix->GetAxis(1)->SetRange(1,fTriggerEventMatrix->GetAxis(1)->GetNbins());\r
+\r
+  //\r
+  // rec. vs true track pt correlation matrix\r
+  //\r
+  hs = (THnSparse*)fTrackPtCorrelationMatrix->Clone("track_pt_correlation_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  //\r
+  // trigger efficiency for INEL\r
+  //\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fGenEventMatrix->Projection(0),"zv_trig_INEL_eff_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+  //\r
+  // trigger bias correction (MB to INEL)\r
+  //\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix,fTriggerEventMatrix,"zv_mult_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0),fTriggerEventMatrix->Projection(0),"zv_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(0,1),fTriggerEventMatrix->Projection(0,1),"zv_mult_trig_MBtoInel_corr_matrix_2D");\r
+  aFolderObj->Add(h2D);\r
+\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fGenEventMatrix->Projection(1),fTriggerEventMatrix->Projection(1),"mult_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+  //\r
+  // event vertex reconstruction correction (MB)\r
+  //\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix,fRecEventMatrix,"zv_mult_event_corr_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0,1),fRecEventMatrix->Projection(0,1),"zv_mult_event_corr_matrix_2D");\r
+  aFolderObj->Add(h2D);\r
+\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(1),fRecEventMatrix->Projection(1),"mult_event_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerEventMatrix->Projection(0),fRecEventMatrix->Projection(0),"zv_event_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+\r
+  //\r
+  // track-event trigger bias correction (MB to INEL)\r
+  //\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix,fTriggerTrackEventMatrix,"zv_pt_eta_track_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,2),fTriggerTrackEventMatrix->Projection(1,2),"eta_pt_track_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(1,0),fTriggerTrackEventMatrix->Projection(1,0),"pt_zv_track_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenTrackEventMatrix->Projection(2,0),fTriggerTrackEventMatrix->Projection(2,0),"zv_eta_track_trig_MBtoInel_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  // efficiency\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1),fGenTrackEventMatrix->Projection(1),"pt_track_trig_MBtoInel_eff_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+  //\r
+  // track-event vertex reconstruction correction (MB)\r
+  //\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix,fRecTrackEventMatrix,"zv_pt_eta_track_event_corr_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,2),fRecTrackEventMatrix->Projection(1,2),"eta_pt_track_event_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(1,0),fRecTrackEventMatrix->Projection(1,0),"pt_zv_track_event_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fTriggerTrackEventMatrix->Projection(2,0),fRecTrackEventMatrix->Projection(2,0),"zv_eta_track_event_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+  \r
+  // efficiency\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecTrackEventMatrix->Projection(1),fTriggerTrackEventMatrix->Projection(1),"pt_track_event_eff_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+  //\r
+  // track rec. efficiency correction\r
+  //\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix,fRecPrimTrackMatrix,"zv_pt_eta_track_corr_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,2),fRecPrimTrackMatrix->Projection(1,2),"eta_pt_track_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1,0),fRecPrimTrackMatrix->Projection(1,0),"pt_zv_track_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2,0),fRecPrimTrackMatrix->Projection(2,0),"zv_eta_track_corr_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  \r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(0),fRecPrimTrackMatrix->Projection(0),"zv_track_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(1),fRecPrimTrackMatrix->Projection(1),"pt_track_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  // efficiency\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecPrimTrackMatrix->Projection(1), fGenPrimTrackMatrix->Projection(1),"pt_track_eff_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fGenPrimTrackMatrix->Projection(2),fRecPrimTrackMatrix->Projection(2),"eta_track_corr_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // secondary track contamination correction\r
+  //\r
+  //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix,fRecTrackMatrix,"zv_pt_eta_track_cont_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_track_cont_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_track_cont_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_track_cont_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_track_cont_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_track_cont_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecSecTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_track_cont_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // multiple track reconstruction correction\r
+  //\r
+  //hs = AlidNdPtHelper::GenerateContCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");\r
+  hs = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix,fRecTrackMatrix,"zv_pt_eta_mult_track_cont_matrix");\r
+  aFolderObj->Add(hs);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,2),fRecTrackMatrix->Projection(1,2),"eta_pt_mult_track_cont_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1,0),fRecTrackMatrix->Projection(1,0),"pt_zv_mult_track_cont_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h2D = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2,0),fRecTrackMatrix->Projection(2,0),"zv_eta_mult_track_cont_matrix");\r
+  aFolderObj->Add(h2D);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(0),fRecTrackMatrix->Projection(0),"zv_mult_track_cont_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(1),fRecTrackMatrix->Projection(1),"pt_mult_track_cont_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  h = AlidNdPtHelper::GenerateCorrMatrix(fRecMultTrackMatrix->Projection(2),fRecTrackMatrix->Projection(2),"eta_mult_track_cont_matrix");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  // Control histograms\r
+  //\r
+  \r
+  if(fHistogramsOn) {\r
+\r
+  // Efficiency electrons, muons, pions, kaons, protons, all\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,1); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,1); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_electrons");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(2,2); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(2,2); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_muons");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(3,3); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(3,3); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_pions");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(4,4); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(4,4); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_kaons");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(5,5); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(5,5); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_protons");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,5); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,5); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_selected");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(2)->SetRange(1,6); \r
+  fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(1,6); \r
+  h1 = fMCPrimTrackHist1[1]->Projection(0);\r
+  h2 = fMCPrimTrackHist1[2]->Projection(0);\r
+  h2c = (TH1D *)h2->Clone();\r
+  h2c->Divide(h1);\r
+  h2c->SetName("eff_pt_all");\r
+  aFolderObj->Add(h2c);\r
+\r
+  fMCPrimTrackHist1[1]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[1]->GetAxis(1)->GetNbins()); \r
+  fMCPrimTrackHist1[2]->GetAxis(1)->SetRange(1,fMCPrimTrackHist1[2]->GetAxis(1)->GetNbins());\r
+\r
+  //  pt spetra\r
+  // - rec, primaries, secondaries\r
+  // - primaries (pid) \r
+  // - secondaries (pid)\r
+  // - secondaries (mech)\r
+  // - secondaries (mother)\r
+  //\r
+\r
+  TH1D *mcPtAccall = fMCTrackHist1[1]->Projection(0);\r
+  mcPtAccall->SetName("mc_pt_acc_all");\r
+  aFolderObj->Add(mcPtAccall);\r
+\r
+  TH1D *mcPtAccprim = fMCPrimTrackHist1[1]->Projection(0);\r
+  mcPtAccprim->SetName("mc_pt_acc_prim");\r
+  aFolderObj->Add(mcPtAccprim);\r
+\r
+  TH1D *mcPtRecall = fMCTrackHist1[2]->Projection(0);\r
+  mcPtRecall->SetName("mc_pt_rec_all");\r
+  aFolderObj->Add(mcPtRecall);\r
+\r
+  TH1D *mcPtRecprim = fMCPrimTrackHist1[2]->Projection(0);\r
+  mcPtRecprim->SetName("mc_pt_rec_prim");\r
+  aFolderObj->Add(mcPtRecprim);\r
+\r
+  TH1D *mcPtRecsec = fMCSecTrackHist1[2]->Projection(0);\r
+  mcPtRecsec->SetName("mc_pt_rec_sec");\r
+  aFolderObj->Add(mcPtRecsec);\r
+\r
+  for(Int_t i = 0; i<6; i++) \r
+  { \r
+    sprintf(name,"mc_pt_rec_prim_pid_%d",i); \r
+    fMCPrimTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);\r
+    h = fMCPrimTrackHist1[2]->Projection(0);\r
+    h->SetName(name);\r
+    aFolderObj->Add(h);\r
+\r
+    sprintf(name,"mc_pt_rec_sec_pid_%d",i); \r
+    fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);\r
+    h = fMCSecTrackHist1[2]->Projection(0);\r
+    h->SetName(name);\r
+    aFolderObj->Add(h);\r
+\r
+    // production mechanisms for given pid\r
+    fMCSecTrackHist1[2]->GetAxis(2)->SetRange(i+1,i+1);\r
+\r
+    for(Int_t j=0; j<20; j++) {\r
+      if(j == 4) {\r
+        // decay\r
+       \r
+        sprintf(name,"mc_pt_rec_sec_pid_%d_decay",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(0);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_eta_rec_sec_pid_%d_decay",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(1);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_mother_rec_sec_pid_%d_decay",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(4);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+      } else if (j == 5) {\r
+       // conversion\r
+\r
+        sprintf(name,"mc_pt_rec_sec_pid_%d_conv",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(0);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_eta_rec_sec_pid_%d_conv",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(1);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_mother_rec_sec_pid_%d_conv",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(4);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+     } else if (j == 13) {\r
+       // mat\r
+       \r
+        sprintf(name,"mc_pt_rec_sec_pid_%d_mat",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(0);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_eta_rec_sec_pid_%d_mat",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(1);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_eta_mother_rec_sec_pid_%d_mat",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(4,1);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_mother_rec_sec_pid_%d_mat",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(4);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+        sprintf(name,"mc_pt_mother_rec_sec_pid_%d_mat",i); \r
+        fMCSecTrackHist1[2]->GetAxis(3)->SetRange(j+1,j+1);\r
+        h = fMCSecTrackHist1[2]->Projection(4,0);\r
+        h->SetName(name);\r
+        aFolderObj->Add(h);\r
+\r
+     } else {\r
+       continue;\r
+     }\r
+   }\r
+\r
+  }\r
+  } // end fHistogramOn\r
+\r
+  //\r
+  //  resolution histograms\r
+  //  only for reconstructed tracks\r
+  //\r
+\r
+  TH2F *h2F=0;\r
+  TCanvas * c = new TCanvas("resol","resol");\r
+  c->cd();\r
+\r
+  //\r
+  fRecMCTrackHist1->GetAxis(1)->SetRangeUser(-0.8,0.79); \r
+\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);\r
+  h->SetXTitle("p_{tmc} (GeV/c)");\r
+  h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");\r
+  h->Draw();\r
+  h->SetName("pt_resolution_vs_mcpt");\r
+  aFolderObj->Add(h);\r
+\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(2,0);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);\r
+  h->SetXTitle("p_{tmc} (GeV/c)");\r
+  h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");\r
+  h->Draw();\r
+  h->SetName("dpt_mean_vs_mcpt");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);\r
+  h->SetXTitle("p_{tmc} (GeV/c)");\r
+  h->SetYTitle("(#eta-#eta_{mc}) resolution");\r
+  h->Draw();\r
+  h->SetName("eta_resolution_vs_mcpt");\r
+  aFolderObj->Add(h);\r
+\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(3,0);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);\r
+  h->SetXTitle("p_{tmc} (GeV/c)");\r
+  h->SetYTitle("(#eta-mc#eta) mean");\r
+  h->Draw();\r
+  h->SetName("deta_mean_vs_mcpt");\r
+  aFolderObj->Add(h);\r
+  \r
+  // \r
+  fRecMCTrackHist1->GetAxis(1)->SetRange(1,fRecMCTrackHist1->GetAxis(1)->GetNbins()); \r
+\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);\r
+  h->SetXTitle("#eta_{mc}");\r
+  h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} resolution");\r
+  h->Draw();\r
+  h->SetName("pt_resolution_vs_mceta");\r
+  aFolderObj->Add(h);\r
+\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(2,1);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);\r
+  h->SetXTitle("#eta_{mc}");\r
+  h->SetYTitle("(p_{t}-p_{tmc})/p_{tmc} mean");\r
+  h->Draw();\r
+  h->SetName("dpt_mean_vs_mceta");\r
+  aFolderObj->Add(h);\r
+\r
+  //\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,0,kTRUE,10);\r
+  h->SetXTitle("#eta_{mc}");\r
+  h->SetYTitle("(#eta-#eta_{mc}) resolution");\r
+  h->Draw();\r
+  h->SetName("eta_resolution_vs_mceta");\r
+  aFolderObj->Add(h);\r
+\r
+  h2F = (TH2F*)fRecMCTrackHist1->Projection(3,1);\r
+  h = AlidNdPtHelper::MakeResol(h2F,1,1,kTRUE,10);\r
+  h->SetXTitle("#eta_{mc}");\r
+  h->SetYTitle("(#eta-mc#eta) mean");\r
+  h->Draw();\r
+  h->SetName("deta_mean_vs_mceta");\r
+  aFolderObj->Add(h);\r
+\r
+  fRecMCTrackHist1->GetAxis(0)->SetRange(1,fRecMCTrackHist1->GetAxis(0)->GetNbins()); \r
+\r
+  } // end use MC info\r
+\r
+  // export objects to analysis folder\r
+  fAnalysisFolder = ExportToFolder(aFolderObj);\r
+\r
+  // delete only TObjArray\r
+  if(aFolderObj) delete aFolderObj;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+TFolder* AlidNdPtAnalysisPbPb::ExportToFolder(TObjArray * const array) \r
+{\r
+  // recreate folder avery time and export objects to new one\r
+  //\r
+  AlidNdPtAnalysisPbPb * comp=this;\r
+  TFolder *folder = comp->GetAnalysisFolder();\r
+\r
+  TString name, title;\r
+  TFolder *newFolder = 0;\r
+  Int_t i = 0;\r
+  Int_t size = array->GetSize();\r
+\r
+  if(folder) { \r
+     // get name and title from old folder\r
+     name = folder->GetName();  \r
+     title = folder->GetTitle();  \r
+\r
+        // delete old one\r
+     delete folder;\r
+\r
+        // create new one\r
+     newFolder = CreateFolder(name.Data(),title.Data());\r
+     newFolder->SetOwner();\r
+\r
+        // add objects to folder\r
+     while(i < size) {\r
+          newFolder->Add(array->At(i));\r
+          i++;\r
+        }\r
+  }\r
+\r
+return newFolder;\r
+}\r
+\r
+//_____________________________________________________________________________\r
+TFolder* AlidNdPtAnalysisPbPb::CreateFolder(TString name,TString title) { \r
+// create folder for analysed histograms\r
+//\r
+TFolder *folder = 0;\r
+  folder = new TFolder(name.Data(),title.Data());\r
+\r
+  return folder;\r
+}\r
diff --git a/PWG0/dNdPt/AlidNdPtAnalysisPbPb.h b/PWG0/dNdPt/AlidNdPtAnalysisPbPb.h
new file mode 100644 (file)
index 0000000..2c29633
--- /dev/null
@@ -0,0 +1,188 @@
+
+#ifndef AlidNdPtAnalysisPbPb_H
+#define AlidNdPtAnalysisPbPb_H
+
+
+//------------------------------------------------------------------------------
+// AlidNdPtAnalysisPbPb class used for dNdPt analysis.in PbPb collision 
+// 
+// Author: J.Otwinowski 04/11/2008 
+//------------------------------------------------------------------------------
+
+class iostream;
+
+class TFile;
+class TCint;
+class TProfile;
+class TFolder;
+class TObjArray;
+class TString;
+class THnSparse;
+
+class AliESDtrackCuts;
+class AliVertexerTracks;
+class AliESD;
+class AliESDfriend;
+class AliESDfriendTrack;
+class AlidNdPtHelper;
+
+#include "AlidNdPt.h"
+
+class AlidNdPtAnalysisPbPb : public AlidNdPt {
+public :
+  AlidNdPtAnalysisPbPb(); 
+  AlidNdPtAnalysisPbPb(Char_t* name, Char_t* title);
+  ~AlidNdPtAnalysisPbPb();
+
+  // Init data members
+  virtual void Init();
+
+  // Process events
+  virtual void Process(AliESDEvent *const esdEvent=0, AliMCEvent *const mcEvent=0);
+
+  // Merge output objects (needed by PROOF) 
+  virtual Long64_t Merge(TCollection* const list);
+
+  // Analyse output histograms 
+  virtual void Analyse();
+
+  // Export objects to folder
+  virtual TFolder *ExportToFolder(TObjArray * const array=0);
+
+  // Get analysis folder
+  TFolder* GetAnalysisFolder() const {return fAnalysisFolder;}
+
+  // Fill control histograms
+  void SetHistogramsOn(const Bool_t histOn=kTRUE) {fHistogramsOn = histOn;}
+  Bool_t IsHistogramsOn() const {return fHistogramsOn;}
+    
+  // Create folder for analysed histograms
+  TFolder *CreateFolder(TString folder = "folderdNdPtAnalysis",TString title = "Analysed dNdPt histograms");
+
+  // Fill histograms
+  void FillHistograms(AliESDtrack *const esdTrack, AliStack *const stack, AlidNdPtHelper::TrackObject trackObj);
+  void FillHistograms(AliStack *const stack, Int_t label, AlidNdPtHelper::TrackObject trackObj);
+  void 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);
+
+  // Getters
+  THnSparseF *GetTrackPtCorrelationMatrix()   const {return fTrackPtCorrelationMatrix;}
+  //
+  THnSparseF *GetGenEventMatrix() const {return fGenEventMatrix;}
+  THnSparseF *GetTriggerEventMatrix() const {return fTriggerEventMatrix;}
+  THnSparseF *GetRecEventMatrix() const {return fRecEventMatrix;}
+  // 
+  THnSparseF *GetGenTrackEventMatrix() const {return fGenTrackEventMatrix;}
+  THnSparseF *GetTriggerTrackEventMatrix() const {return fTriggerTrackEventMatrix;}
+  THnSparseF *GetRecTrackEventMatrix() const {return fRecTrackEventMatrix;}
+  //
+  THnSparseF *GetGenTrackMatrix() const {return fGenTrackMatrix;}
+  THnSparseF *GetGenPrimTrackMatrix() const {return fGenPrimTrackMatrix;}
+  THnSparseF *GetRecPrimTrackMatrix() const {return fRecPrimTrackMatrix;}
+
+  THnSparseF *GetRecTrackMatrix() const {return fRecTrackMatrix;}
+  THnSparseF *GetRecSecTrackMatrix() const {return fRecSecTrackMatrix;}
+  THnSparseF *GetRecMultTrackMatrix() const {return fRecMultTrackMatrix;}
+  //
+  // control histograms
+  //
+  THnSparseF *GetMCEventHist1() const {return fMCEventHist1;}
+  THnSparseF *GetRecEventHist1() const {return fRecEventHist1;}
+  THnSparseF *GetRecEventHist2() const {return fRecEventHist2;}
+
+
+  THnSparseF *GetRecMCEventHist1() const {return fRecMCEventHist1;}
+  THnSparseF *GetRecMCTrackHist1() const {return fRecMCTrackHist1;}
+
+  THnSparseF *GetRecMCEventHist2() const {return fRecMCEventHist2;}
+
+  THnSparseF *GetMCTrackHist1(Int_t i) const {return fMCTrackHist1[i];}
+  THnSparseF *GetMCPrimTrackHist1(Int_t i) const {return fMCPrimTrackHist1[i];}
+  THnSparseF *GetMCPrimTrackHist2(Int_t i) const {return fMCPrimTrackHist2[i];}
+  THnSparseF *GetMCSecTrackHist1(Int_t i) const {return fMCSecTrackHist1[i];}
+
+  THnSparseF *GetRecTrackHist1(Int_t i) const {return fRecTrackHist1[i];}
+  THnSparseF *GetRecTrackMultHist1(Int_t i) const {return fRecTrackMultHist1[i];}
+
+
+  THnSparseF *GetMCMultRecTrackHist1() const {return fMCMultRecTrackHist1;}
+
+  THnSparseF *GetRecTrackHist2() const {return fRecTrackHist2;}
+
+private:
+
+  // analysis folder 
+  TFolder *fAnalysisFolder; // folder for analysed histograms
+  Bool_t fHistogramsOn; // switch on/off filling of control histograms 
+
+  // 
+  // correlation matrices (histograms)
+  //
+  // rec. track pt vs true track pt correlation matrix for given eta
+  THnSparseF *fTrackPtCorrelationMatrix; //-> Pt:mcPt:mcEta
+
+  //
+  // event level correction 
+  //
+  // all generated
+  THnSparseF *fGenEventMatrix; //-> mcZv:multMB 
+
+  // trigger bias corrections (fTriggerEventMatrix / fGenEventMatrix)
+  THnSparseF *fTriggerEventMatrix; //-> mcZv:multMB
+
+  // event vertex rec. eff correction (fRecEventMatrix / fTriggerEventMatrix)
+  THnSparseF *fRecEventMatrix; //-> mcZv:multMB 
+
+  // track-event level correction 
+  THnSparseF *fGenTrackEventMatrix; //-> mcZv:mcPt:mcEta
+
+  // trigger bias corrections (fTriggerTrackEventMatrix / fGenTrackEventMatrix)
+  THnSparseF *fTriggerTrackEventMatrix; //-> mcZv:mcPt:mcEta
+
+  // event vertex rec. corrections (fRecTrackEventMatrix / fTriggerTrackEventMatrix)
+  THnSparseF *fRecTrackEventMatrix; //-> mcZv:Pt:mcEta
+
+  //
+  // track level correction 
+  //
+  // track rec. efficiency correction (fRecPrimTrackMatrix / fGenPrimTrackMatrix)
+  THnSparseF *fGenTrackMatrix; //-> mcZv:mcPt:mcEta
+  THnSparseF *fGenPrimTrackMatrix; //-> mcZv:mcPt:mcEta
+  THnSparseF *fRecPrimTrackMatrix; //-> mcZv:mcPt:mcEta
+  // secondary track contamination correction (fRecSecTrackMatrix / fRecTrackMatrix)
+  THnSparseF *fRecTrackMatrix;    //-> mcZv:mcPt:mcEta
+  THnSparseF *fRecSecTrackMatrix; //-> mcZv:mcPt:mcEta
+  // multiple rec. track corrections (fRecMultTrackMatrix / fRecTrackMatrix)
+  THnSparseF *fRecMultTrackMatrix; //-> mcZv:Pt:mcEta
+
+  //
+  // ESD and MC control analysis histograms
+  //
+  // THnSparse event histograms
+  THnSparseF *fMCEventHist1;  //-> mcXv:mcYv:mcZv
+  THnSparseF *fRecEventHist1; //-> Xv:Yv:Zv
+  THnSparseF *fRecEventHist2; //-> Zv:multMB:mult
+  THnSparseF *fRecMCEventHist1; //-> Xv-mcXv:Yv-mcYv:Zv-mcZv
+  THnSparseF *fRecMCEventHist2; //-> Xv-mcXv:Zv-mcZv:mult
+
+  // [0] - after charged track selection, [1] - after acceptance cuts, [2] - after esd track cuts
+  THnSparseF *fMCTrackHist1[AlidNdPtHelper::kCutSteps];     //-> mcPt:mcEta:mcPhi
+  THnSparseF *fMCPrimTrackHist1[AlidNdPtHelper::kCutSteps]; //-> mcPt:mcEta:pid:mech:mother
+  THnSparseF *fMCPrimTrackHist2[AlidNdPtHelper::kCutSteps]; //-> pdg:mech:mother
+  THnSparseF *fMCSecTrackHist1[AlidNdPtHelper::kCutSteps];  //-> mcPt:mcEta:pid:mech:mother
+
+  THnSparseF *fRecTrackHist1[AlidNdPtHelper::kCutSteps];     //-> Pt:Eta:Phi
+  THnSparseF *fRecTrackMultHist1[AlidNdPtHelper::kCutSteps]; //-> Pt:mult
+  THnSparseF *fRecMCTrackHist1; //-> mcPt:mcEta:(Pt-mcPt)/mcPt:(Eta-mcEta)
+
+  //multple reconstructed tracks
+  THnSparseF *fMCMultRecTrackHist1; //-> mcPt:mcEta:pid
+  // track control histograms
+  THnSparseF *fRecTrackHist2;  //-> nclust:chi2:Pt:Eta:Phi
+
+  AlidNdPtAnalysisPbPb(const AlidNdPtAnalysisPbPb&); // not implemented
+  AlidNdPtAnalysisPbPb& operator=(const AlidNdPtAnalysisPbPb&); // not implemented
+
+  ClassDef(AlidNdPtAnalysisPbPb,2);
+};
+
+#endif
index 3897d12..67432b7 100644 (file)
@@ -6,6 +6,7 @@
 
 SRCS = dNdPt/AlidNdPtHelper.cxx \
       dNdPt/AlidNdPtAnalysis.cxx \
+      dNdPt/AlidNdPtAnalysisPbPb.cxx \
       dNdPt/AlidNdPtCorrection.cxx \
       dNdPt/AlidNdPtAcceptanceCuts.cxx \
       dNdPt/AlidNdPtEventCuts.cxx \