]> git.uio.no Git - u/mrichter/AliRoot.git/blob - PWGLF/SPECTRA/ChargedHadrons/dNdPt/AlidNdPtCorrection.cxx
changes for Vertex and Tracks classes
[u/mrichter/AliRoot.git] / PWGLF / SPECTRA / ChargedHadrons / dNdPt / AlidNdPtCorrection.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 // AlidNdPtCorrection class:
17 //
18 // a. functionality:
19 // - applies corrections on dNdPt spectra
20 // - fills corrected dNdPt histograms
21 // - fills correction control histograms 
22 //
23 // b. data members:
24 // - dNdPt spectra before and after correction procedure
25 // - control histograms
26 // - correction matrices (must be loaded)
27 // 
28 // Author: J.Otwinowski 04/11/2008 
29 //------------------------------------------------------------------------------
30
31 #include "TFile.h"
32 #include "TH1.h"
33 #include "TH2.h"
34
35 #include "AliHeader.h"  
36 #include "AliGenEventHeader.h"  
37 #include "AliStack.h"  
38 #include "AliESDEvent.h"  
39 #include "AliMCEvent.h"  
40 #include "AliESDtrackCuts.h"  
41 #include "AliLog.h" 
42 #include "AliMultiplicity.h"
43 #include "AliTracker.h"
44
45 #include "AlidNdPtEventCuts.h"
46 #include "AlidNdPtAcceptanceCuts.h"
47 #include "AliPhysicsSelection.h"
48
49 #include "AliPWG0Helper.h"
50 #include "AlidNdPtHelper.h"
51 #include "AlidNdPtAnalysis.h"
52 #include "AlidNdPtCorrection.h"
53
54 using namespace std;
55
56 ClassImp(AlidNdPtCorrection)
57
58 //_____________________________________________________________________________
59 //AlidNdPtCorrection::AlidNdPtCorrection(): TNamed(),
60   AlidNdPtCorrection::AlidNdPtCorrection(): AlidNdPt(),
61   fCorrectionFolder(0),
62   fMCEventHist1(0),
63   fRecEventHist1(0),
64   fRecEventMultHist1(0),
65   fMCAllEventMultHist1(0),
66   fMCAllNDEventMultHist1(0),
67   fMCAllNSDEventMultHist1(0),
68   fMCTriggerMultHist1(0),
69   fMCEventMultHist1(0),
70   fMCAllPrimTrackMultHist1(0),
71   fMCNDEventAllPrimTrackMultHist1(0),
72   fMCNSDEventAllPrimTrackMultHist1(0),
73   fMCTriggerPrimTrackMultHist1(0),
74   fMCEventPrimTrackMultHist1(0),
75   fMCAllPrimTrackTrueMultHist1(0),
76   fMCNDEventAllPrimTrackTrueMultHist1(0),
77   fMCNSDEventAllPrimTrackTrueMultHist1(0),
78   fMCTriggerPrimTrackTrueMultHist1(0),
79   fMCEventPrimTrackTrueMultHist1(0),
80   fMCAllPrimTrackTrueMultHist2(0),
81   fMCNDEventAllPrimTrackTrueMultHist2(0),
82   fMCNSDEventAllPrimTrackTrueMultHist2(0),
83   fMCTriggerPrimTrackTrueMultHist2(0),
84   fMCEventPrimTrackTrueMultHist2(0),
85   fMCAllPrimTrackMeanPtMult1(0),
86   fMCNDEventAllPrimTrackMeanPtMult1(0),
87   fMCNSDEventAllPrimTrackMeanPtMult1(0),
88   fMCTriggerPrimTrackMeanPtMult1(0),
89   fMCEventPrimTrackMeanPtMult1(0),
90   fMCAllPrimTrackMeanPtTrueMult1(0),
91   fMCNDEventAllPrimTrackMeanPtTrueMult1(0),
92   fMCNSDEventAllPrimTrackMeanPtTrueMult1(0),
93   fMCTriggerPrimTrackMeanPtTrueMult1(0),
94   fMCEventPrimTrackMeanPtTrueMult1(0),
95   fEventMultCorrelationMatrix(0),
96   fZvNorm(0),
97   fZvEmptyEventsNorm(0),
98   fLHCBin0Background(0),
99   fCorrTriggerMBtoInelEventMatrix(0),
100   fCorrTriggerMBtoNDEventMatrix(0),
101   fCorrTriggerMBtoNSDEventMatrix(0),
102   fCorrEventMatrix(0),
103   fCorrTriggerMBtoInelTrackEventMatrix(0),
104   fCorrTriggerMBtoNDTrackEventMatrix(0),
105   fCorrTriggerMBtoNSDTrackEventMatrix(0),
106   fCorrTrackEventMatrix(0),
107   fCorrTrackMatrix(0),
108   fCorrHighPtTrackMatrix(0),
109   fContTrackMatrix(0),
110   fContMultTrackMatrix(0),
111   fCorrMatrixFileName(""),
112   fCosmicsHisto(0),
113   fEventCount(0)
114 {
115   // default constructor
116   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { 
117     fRecTrackHist1[i]=0;     
118   }
119
120   for(Int_t i=0; i<8; i++) { 
121     fCorrRecTrackMultHist1[i] = 0;
122     fCorrRecTrackTrueMultHist1[i] = 0;
123     fCorrRecTrackTrueMultHist2[i] = 0;
124     fCorrRecTrackMeanPtMultHist1[i] = 0;
125     fCorrRecTrackMeanPtTrueMultHist1[i] = 0;
126     fCorrRecTrackPt1[i] = 0;
127   }
128
129   for(Int_t i=0; i<5; i++) { 
130     fCorrRecEventHist1[i] = 0;
131     fCorrRecEventHist2[i] = 0;
132   }
133
134   Init();
135 }
136
137 //_____________________________________________________________________________
138 AlidNdPtCorrection::AlidNdPtCorrection(Char_t* name, Char_t* title, TString corrMatrixFileName): AlidNdPt(name,title),
139   fCorrectionFolder(0),
140   fMCEventHist1(0),
141   fRecEventHist1(0),
142   fRecEventMultHist1(0),
143   fMCAllEventMultHist1(0),
144   fMCAllNDEventMultHist1(0),
145   fMCAllNSDEventMultHist1(0),
146   fMCTriggerMultHist1(0),
147   fMCEventMultHist1(0),
148   fMCAllPrimTrackMultHist1(0),
149   fMCNDEventAllPrimTrackMultHist1(0),
150   fMCNSDEventAllPrimTrackMultHist1(0),
151   fMCTriggerPrimTrackMultHist1(0),
152   fMCEventPrimTrackMultHist1(0),
153   fMCAllPrimTrackTrueMultHist1(0),
154   fMCNDEventAllPrimTrackTrueMultHist1(0),
155   fMCNSDEventAllPrimTrackTrueMultHist1(0),
156   fMCTriggerPrimTrackTrueMultHist1(0),
157   fMCEventPrimTrackTrueMultHist1(0),
158   fMCAllPrimTrackTrueMultHist2(0),
159   fMCNDEventAllPrimTrackTrueMultHist2(0),
160   fMCNSDEventAllPrimTrackTrueMultHist2(0),
161   fMCTriggerPrimTrackTrueMultHist2(0),
162   fMCEventPrimTrackTrueMultHist2(0),
163   fMCAllPrimTrackMeanPtMult1(0),
164   fMCNDEventAllPrimTrackMeanPtMult1(0),
165   fMCNSDEventAllPrimTrackMeanPtMult1(0),
166   fMCTriggerPrimTrackMeanPtMult1(0),
167   fMCEventPrimTrackMeanPtMult1(0),
168   fMCAllPrimTrackMeanPtTrueMult1(0),
169   fMCNDEventAllPrimTrackMeanPtTrueMult1(0),
170   fMCNSDEventAllPrimTrackMeanPtTrueMult1(0),
171   fMCTriggerPrimTrackMeanPtTrueMult1(0),
172   fMCEventPrimTrackMeanPtTrueMult1(0),
173   fEventMultCorrelationMatrix(0),
174   fZvNorm(0),
175   fZvEmptyEventsNorm(0),
176   fLHCBin0Background(0),
177   fCorrTriggerMBtoInelEventMatrix(0),
178   fCorrTriggerMBtoNDEventMatrix(0),
179   fCorrTriggerMBtoNSDEventMatrix(0),
180   fCorrEventMatrix(0),
181   fCorrTriggerMBtoInelTrackEventMatrix(0),
182   fCorrTriggerMBtoNDTrackEventMatrix(0),
183   fCorrTriggerMBtoNSDTrackEventMatrix(0),
184   fCorrTrackEventMatrix(0),
185   fCorrTrackMatrix(0),
186   fCorrHighPtTrackMatrix(0),
187   fContTrackMatrix(0),
188   fContMultTrackMatrix(0),
189   fCorrMatrixFileName(corrMatrixFileName),
190   fCosmicsHisto(0),
191   fEventCount(0)
192 {
193   // constructor
194   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { 
195     fRecTrackHist1[i]=0;     
196   }
197
198   for(Int_t i=0; i<8; i++) { 
199     fCorrRecTrackMultHist1[i] = 0;
200     fCorrRecTrackTrueMultHist1[i] = 0;
201     fCorrRecTrackTrueMultHist2[i] = 0;
202     fCorrRecTrackMeanPtMultHist1[i] = 0;
203     fCorrRecTrackMeanPtTrueMultHist1[i] = 0;
204     fCorrRecTrackPt1[i] = 0;
205   }
206
207   for(Int_t i=0; i<5; i++) { 
208     fCorrRecEventHist1[i] = 0;
209     fCorrRecEventHist2[i] = 0;
210   }
211
212   Init();
213 }
214
215 //_____________________________________________________________________________
216 AlidNdPtCorrection::~AlidNdPtCorrection() {
217   // 
218   // destructor
219   //
220   if(fMCEventHist1) delete fMCEventHist1; fMCEventHist1=0;
221   if(fRecEventHist1) delete fRecEventHist1; fRecEventHist1=0;
222   if(fRecEventMultHist1) delete fRecEventMultHist1; fRecEventMultHist1=0;
223
224   if(fMCAllEventMultHist1) delete fMCAllEventMultHist1; fMCAllEventMultHist1=0;
225   if(fMCAllNDEventMultHist1) delete fMCAllNDEventMultHist1; fMCAllNDEventMultHist1=0;
226   if(fMCAllNSDEventMultHist1) delete fMCAllNSDEventMultHist1; fMCAllNSDEventMultHist1=0;
227   if(fMCTriggerMultHist1) delete fMCTriggerMultHist1; fMCTriggerMultHist1=0;
228   if(fMCEventMultHist1) delete fMCEventMultHist1; fMCEventMultHist1=0;
229
230   if(fMCAllPrimTrackMultHist1) delete fMCAllPrimTrackMultHist1; fMCAllPrimTrackMultHist1=0;
231   if(fMCNDEventAllPrimTrackMultHist1) delete fMCNDEventAllPrimTrackMultHist1; fMCNDEventAllPrimTrackMultHist1=0;
232   if(fMCNSDEventAllPrimTrackMultHist1) delete fMCNSDEventAllPrimTrackMultHist1; fMCNSDEventAllPrimTrackMultHist1=0;
233   if(fMCTriggerPrimTrackMultHist1) delete fMCTriggerPrimTrackMultHist1; fMCTriggerPrimTrackMultHist1=0;
234   if(fMCEventPrimTrackMultHist1) delete fMCEventPrimTrackMultHist1; fMCEventPrimTrackMultHist1=0;
235
236   if(fMCAllPrimTrackTrueMultHist1) delete fMCAllPrimTrackTrueMultHist1; fMCAllPrimTrackTrueMultHist1=0;
237   if(fMCNDEventAllPrimTrackTrueMultHist1) delete fMCNDEventAllPrimTrackTrueMultHist1; fMCNDEventAllPrimTrackTrueMultHist1=0;
238   if(fMCNSDEventAllPrimTrackTrueMultHist1) delete fMCNSDEventAllPrimTrackTrueMultHist1; fMCNSDEventAllPrimTrackTrueMultHist1=0;
239   if(fMCTriggerPrimTrackTrueMultHist1) delete fMCTriggerPrimTrackTrueMultHist1; fMCTriggerPrimTrackTrueMultHist1=0;
240   if(fMCEventPrimTrackTrueMultHist1) delete fMCEventPrimTrackTrueMultHist1; fMCEventPrimTrackTrueMultHist1=0;
241
242   if(fMCAllPrimTrackTrueMultHist2) delete fMCAllPrimTrackTrueMultHist2; fMCAllPrimTrackTrueMultHist2=0;
243   if(fMCNDEventAllPrimTrackTrueMultHist2) delete fMCNDEventAllPrimTrackTrueMultHist2; fMCNDEventAllPrimTrackTrueMultHist2=0;
244   if(fMCNSDEventAllPrimTrackTrueMultHist2) delete fMCNSDEventAllPrimTrackTrueMultHist2; fMCNSDEventAllPrimTrackTrueMultHist2=0;
245   if(fMCTriggerPrimTrackTrueMultHist2) delete fMCTriggerPrimTrackTrueMultHist2; fMCTriggerPrimTrackTrueMultHist2=0;
246   if(fMCEventPrimTrackTrueMultHist2) delete fMCEventPrimTrackTrueMultHist2; fMCEventPrimTrackTrueMultHist2=0;
247
248
249
250
251   if(fMCAllPrimTrackMeanPtMult1) delete fMCAllPrimTrackMeanPtMult1; fMCAllPrimTrackMeanPtMult1=0;
252   if(fMCNDEventAllPrimTrackMeanPtMult1) delete fMCNDEventAllPrimTrackMeanPtMult1; fMCNDEventAllPrimTrackMeanPtMult1=0;
253   if(fMCNSDEventAllPrimTrackMeanPtMult1) delete fMCNSDEventAllPrimTrackMeanPtMult1; fMCNSDEventAllPrimTrackMeanPtMult1=0;
254   if(fMCTriggerPrimTrackMeanPtMult1) delete fMCTriggerPrimTrackMeanPtMult1; fMCTriggerPrimTrackMeanPtMult1=0;
255   if(fMCEventPrimTrackMeanPtMult1) delete fMCEventPrimTrackMeanPtMult1; fMCEventPrimTrackMeanPtMult1=0;
256
257   if(fMCAllPrimTrackMeanPtTrueMult1) delete fMCAllPrimTrackMeanPtTrueMult1; fMCAllPrimTrackMeanPtTrueMult1=0;
258   if(fMCNDEventAllPrimTrackMeanPtTrueMult1) delete fMCNDEventAllPrimTrackMeanPtTrueMult1; fMCNDEventAllPrimTrackMeanPtTrueMult1=0;
259   if(fMCNSDEventAllPrimTrackMeanPtTrueMult1) delete fMCNSDEventAllPrimTrackMeanPtTrueMult1; fMCNSDEventAllPrimTrackMeanPtTrueMult1=0;
260   if(fMCTriggerPrimTrackMeanPtTrueMult1) delete fMCTriggerPrimTrackMeanPtTrueMult1; fMCTriggerPrimTrackMeanPtTrueMult1=0;
261   if(fMCEventPrimTrackMeanPtTrueMult1) delete fMCEventPrimTrackMeanPtTrueMult1; fMCEventPrimTrackMeanPtTrueMult1=0;
262
263   if(fCosmicsHisto) delete fCosmicsHisto; fCosmicsHisto=0;
264   if(fEventCount) delete fEventCount; fEventCount=0;
265
266   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { 
267     if(fRecTrackHist1[i]) delete fRecTrackHist1[i]; fRecTrackHist1[i]=0;
268   }
269
270   for(Int_t i=0; i<8; i++) { 
271     if(fCorrRecTrackMultHist1[i]) delete fCorrRecTrackMultHist1[i]; fCorrRecTrackMultHist1[i]=0;
272     if(fCorrRecTrackTrueMultHist1[i]) delete fCorrRecTrackTrueMultHist1[i]; fCorrRecTrackTrueMultHist1[i]=0;
273     if(fCorrRecTrackTrueMultHist2[i]) delete fCorrRecTrackTrueMultHist2[i]; fCorrRecTrackTrueMultHist2[i]=0;
274     if(fCorrRecTrackMeanPtMultHist1[i]) delete fCorrRecTrackMeanPtMultHist1[i]; fCorrRecTrackMeanPtMultHist1[i]=0;
275     if(fCorrRecTrackMeanPtTrueMultHist1[i]) delete fCorrRecTrackMeanPtTrueMultHist1[i]; fCorrRecTrackMeanPtTrueMultHist1[i]=0;
276     if(fCorrRecTrackPt1[i]) delete fCorrRecTrackPt1[i]; fCorrRecTrackPt1[i]=0;
277   }
278
279   for(Int_t i=0; i<5; i++) { 
280     if(fCorrRecEventHist1[i]) delete fCorrRecEventHist1[i]; fCorrRecEventHist1[i]=0;
281     if(fCorrRecEventHist2[i]) delete fCorrRecEventHist2[i]; fCorrRecEventHist2[i]=0;
282   }
283
284   if(fCorrectionFolder) delete fCorrectionFolder; fCorrectionFolder=0;
285 }
286
287 //_____________________________________________________________________________
288 void AlidNdPtCorrection::Init(){
289   //
290   // Init histograms
291   //
292   const Int_t etaNbins = 30; 
293   const Int_t zvNbins = 12;
294
295   // UA1 bining
296   //const Int_t ptNbins = 52; 
297   //Double_t binsPt[ptNbins+1] = { 0.20, 0.30, 0.40, 0.50, 0.60, 0.70, 0.80, 0.90, 1.00, 1.10, 1.20, 1.30, 1.40, 1.50, 1.60, 1.70, 1.80, 1.90, 2.00, 2.10, 2.20, 2.30, 2.40, 2.50, 2.60, 2.70, 2.80, 2.90, 3.00, 3.10, 3.20, 3.30, 3.40, 3.50, 3.60, 3.70, 3.80, 3.90, 4.00, 4.20, 4.40, 4.60, 4.80, 5.00, 5.20, 5.40, 5.60, 5.80, 6.00, 7.00, 7.60, 8.80, 9.60 }; 
298
299  /*
300   const Int_t ptNbins = 56; 
301   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,7.5,8.0,9.0,10.0,11.0,12.0,13.0,14.0,15.0,16.0};
302   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};
303   Double_t binsZv[zvNbins+1] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
304   */
305
306
307   const Int_t ptNbins = 68;
308   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,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};
309
310
311   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};
312   Double_t binsZv[zvNbins+1] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};
313
314
315   //
316   Int_t binsMCEventHist1[3]={100,100,140};
317   Double_t minMCEventHist1[3]={-0.1,-0.1,-35.}; 
318   Double_t maxMCEventHist1[3]={0.1,0.1,35.}; 
319   fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv",3,binsMCEventHist1,minMCEventHist1,maxMCEventHist1);
320   fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)");
321   fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)");
322   fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)");
323   fMCEventHist1->Sumw2();
324
325   //
326   Int_t binsRecEventHist1[3]={100,100,140};
327   Double_t minRecEventHist1[3]={-3.,-3.,-35.}; 
328   Double_t maxRecEventHist1[3]={3.,3.,35.}; 
329   fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv",3,binsRecEventHist1,minRecEventHist1,maxRecEventHist1);
330   fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)");
331   fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)");
332   fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)");
333   fRecEventHist1->Sumw2();
334
335   //
336   Int_t binsRecEventMultHist1[2]={150,150};
337   Double_t minRecEventMultHist1[2]={-0.5,-0.5}; 
338   Double_t maxRecEventMultHist1[2]={149.5,149.5}; 
339   fRecEventMultHist1 = new THnSparseF("fRecEventMultHist1","track multiplicity:tracklet multiplicity",2,binsRecEventMultHist1,minRecEventMultHist1,maxRecEventMultHist1);
340   fRecEventMultHist1->GetAxis(0)->SetTitle("track_mult");
341   fRecEventMultHist1->GetAxis(1)->SetTitle("tracklet_mult");
342   fRecEventMultHist1->Sumw2();
343
344   //
345   char name[256];
346   char title[256];
347
348   Int_t binsMCAllPrimTrackMultHist1[3]={ptNbins,etaNbins,150};
349   Double_t minMCAllPrimTrackMultHist1[3]={0.,-1.,-0.5}; 
350   Double_t maxMCAllPrimTrackMultHist1[3]={20.,1.,149.5}; 
351   snprintf(name,256,"fMCAllPrimTrackMultHist1");
352   snprintf(title,256,"mcPt:mcEta:multiplicity");
353   
354   fMCAllPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCAllPrimTrackMultHist1,minMCAllPrimTrackMultHist1,maxMCAllPrimTrackMultHist1);
355   fMCAllPrimTrackMultHist1->SetBinEdges(0,binsPt);
356   fMCAllPrimTrackMultHist1->SetBinEdges(1,binsEta);
357   fMCAllPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
358   fMCAllPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");
359   fMCAllPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");
360   fMCAllPrimTrackMultHist1->Sumw2();
361
362   Int_t binsMCNDEventAllPrimTrackMultHist1[3]={ptNbins,etaNbins,150};
363   Double_t minMCNDEventAllPrimTrackMultHist1[3]={0.,-1.,-0.5}; 
364   Double_t maxMCNDEventAllPrimTrackMultHist1[3]={20.,1.,149.5}; 
365   snprintf(name,256,"fMCNDEventAllPrimTrackMultHist1");
366   snprintf(title,256,"mcPt:mcEta:multiplicity");
367   
368   fMCNDEventAllPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCNDEventAllPrimTrackMultHist1,minMCNDEventAllPrimTrackMultHist1,maxMCNDEventAllPrimTrackMultHist1);
369   fMCNDEventAllPrimTrackMultHist1->SetBinEdges(0,binsPt);
370   fMCNDEventAllPrimTrackMultHist1->SetBinEdges(1,binsEta);
371   fMCNDEventAllPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
372   fMCNDEventAllPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");
373   fMCNDEventAllPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");
374   fMCNDEventAllPrimTrackMultHist1->Sumw2();
375
376   Int_t binsMCNSDEventAllPrimTrackMultHist1[3]={ptNbins,etaNbins,150};
377   Double_t minMCNSDEventAllPrimTrackMultHist1[3]={0.,-1.,-0.5}; 
378   Double_t maxMCNSDEventAllPrimTrackMultHist1[3]={20.,1.,149.5}; 
379   snprintf(name,256,"fMCNSDEventAllPrimTrackMultHist1");
380   snprintf(title,256,"mcPt:mcEta:multiplicity");
381   
382   fMCNSDEventAllPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCNSDEventAllPrimTrackMultHist1,minMCNSDEventAllPrimTrackMultHist1,maxMCNSDEventAllPrimTrackMultHist1);
383   fMCNSDEventAllPrimTrackMultHist1->SetBinEdges(0,binsPt);
384   fMCNSDEventAllPrimTrackMultHist1->SetBinEdges(1,binsEta);
385   fMCNSDEventAllPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
386   fMCNSDEventAllPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");
387   fMCNSDEventAllPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");
388   fMCNSDEventAllPrimTrackMultHist1->Sumw2();
389
390   Int_t binsMCEventTriggerPrimTrackMultHist1[3]={ptNbins,etaNbins,150};
391   Double_t minMCEventTriggerPrimTrackMultHist1[3]={0.,-1.,-0.5}; 
392   Double_t maxMCEventTriggerPrimTrackMultHist1[3]={20.,1.,149.5}; 
393   snprintf(name,256,"fMCTriggerPrimTrackMultHist1");
394   snprintf(title,256,"mcPt:mcEta:multiplicity");
395   
396   fMCTriggerPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCEventTriggerPrimTrackMultHist1,minMCEventTriggerPrimTrackMultHist1,maxMCEventTriggerPrimTrackMultHist1);
397   fMCTriggerPrimTrackMultHist1->SetBinEdges(0,binsPt);
398   fMCTriggerPrimTrackMultHist1->SetBinEdges(1,binsEta);
399   fMCTriggerPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
400   fMCTriggerPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");
401   fMCTriggerPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");
402   fMCTriggerPrimTrackMultHist1->Sumw2();
403
404   Int_t binsMCEventPrimTrackMultHist1[3]={ptNbins,etaNbins,150};
405   Double_t minMCEventPrimTrackMultHist1[3]={0.,-1.,-0.5}; 
406   Double_t maxMCEventPrimTrackMultHist1[3]={20.,1.,149.5}; 
407   snprintf(name,256,"fMCEventPrimTrackMultHist1");
408   snprintf(title,256,"mcPt:mcEta:multiplicity");
409   
410   fMCEventPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCEventPrimTrackMultHist1,minMCEventPrimTrackMultHist1,maxMCEventPrimTrackMultHist1);
411   fMCEventPrimTrackMultHist1->SetBinEdges(0,binsPt);
412   fMCEventPrimTrackMultHist1->SetBinEdges(1,binsEta);
413   fMCEventPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
414   fMCEventPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");
415   fMCEventPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");
416   fMCEventPrimTrackMultHist1->Sumw2();
417
418   //
419   // true multiplicity
420   //
421
422   Int_t binsMCAllPrimTrackTrueMultHist1[3]={ptNbins,etaNbins,150};
423   Double_t minMCAllPrimTrackTrueMultHist1[3]={0.,-1.,-0.5}; 
424   Double_t maxMCAllPrimTrackTrueMultHist1[3]={20.,1.,149.5}; 
425   snprintf(name,256,"fMCAllPrimTrackTrueMultHist1");
426   snprintf(title,256,"mcPt:mcEta:true_mult");
427   
428   fMCAllPrimTrackTrueMultHist1 = new THnSparseF(name,title,3,binsMCAllPrimTrackTrueMultHist1,minMCAllPrimTrackTrueMultHist1,maxMCAllPrimTrackTrueMultHist1);
429   fMCAllPrimTrackTrueMultHist1->SetBinEdges(0,binsPt);
430   fMCAllPrimTrackTrueMultHist1->SetBinEdges(1,binsEta);
431   fMCAllPrimTrackTrueMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
432   fMCAllPrimTrackTrueMultHist1->GetAxis(1)->SetTitle("mcEta");
433   fMCAllPrimTrackTrueMultHist1->GetAxis(2)->SetTitle("true_mult");
434   fMCAllPrimTrackTrueMultHist1->Sumw2();
435
436   Int_t binsMCNDEventAllPrimTrackTrueMultHist1[3]={ptNbins,etaNbins,150};
437   Double_t minMCNDEventAllPrimTrackTrueMultHist1[3]={0.,-1.,-0.5}; 
438   Double_t maxMCNDEventAllPrimTrackTrueMultHist1[3]={20.,1.,149.5}; 
439   snprintf(name,256,"fMCNDEventAllPrimTrackTrueMultHist1");
440   snprintf(title,256,"mcPt:mcEta:true_mult");
441   
442   fMCNDEventAllPrimTrackTrueMultHist1 = new THnSparseF(name,title,3,binsMCNDEventAllPrimTrackTrueMultHist1,minMCNDEventAllPrimTrackTrueMultHist1,maxMCNDEventAllPrimTrackTrueMultHist1);
443   fMCNDEventAllPrimTrackTrueMultHist1->SetBinEdges(0,binsPt);
444   fMCNDEventAllPrimTrackTrueMultHist1->SetBinEdges(1,binsEta);
445   fMCNDEventAllPrimTrackTrueMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
446   fMCNDEventAllPrimTrackTrueMultHist1->GetAxis(1)->SetTitle("mcEta");
447   fMCNDEventAllPrimTrackTrueMultHist1->GetAxis(2)->SetTitle("true_mult");
448   fMCNDEventAllPrimTrackTrueMultHist1->Sumw2();
449
450   Int_t binsMCNSDEventAllPrimTrackTrueMultHist1[3]={ptNbins,etaNbins,150};
451   Double_t minMCNSDEventAllPrimTrackTrueMultHist1[3]={0.,-1.,-0.5}; 
452   Double_t maxMCNSDEventAllPrimTrackTrueMultHist1[3]={20.,1.,149.5}; 
453   snprintf(name,256,"fMCNSDEventAllPrimTrackTrueMultHist1");
454   snprintf(title,256,"mcPt:mcEta:true_mult");
455   
456   fMCNSDEventAllPrimTrackTrueMultHist1 = new THnSparseF(name,title,3,binsMCNSDEventAllPrimTrackTrueMultHist1,minMCNSDEventAllPrimTrackTrueMultHist1,maxMCNSDEventAllPrimTrackTrueMultHist1);
457   fMCNSDEventAllPrimTrackTrueMultHist1->SetBinEdges(0,binsPt);
458   fMCNSDEventAllPrimTrackTrueMultHist1->SetBinEdges(1,binsEta);
459   fMCNSDEventAllPrimTrackTrueMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
460   fMCNSDEventAllPrimTrackTrueMultHist1->GetAxis(1)->SetTitle("mcEta");
461   fMCNSDEventAllPrimTrackTrueMultHist1->GetAxis(2)->SetTitle("true_mult");
462   fMCNSDEventAllPrimTrackTrueMultHist1->Sumw2();
463
464   Int_t binsMCEventTriggerPrimTrackTrueMultHist1[3]={ptNbins,etaNbins,150};
465   Double_t minMCEventTriggerPrimTrackTrueMultHist1[3]={0.,-1.,-0.5}; 
466   Double_t maxMCEventTriggerPrimTrackTrueMultHist1[3]={20.,1.,149.5}; 
467   snprintf(name,256,"fMCTriggerPrimTrackTrueMultHist1");
468   snprintf(title,256,"mcPt:mcEta:true_mult");
469   
470   fMCTriggerPrimTrackTrueMultHist1 = new THnSparseF(name,title,3,binsMCEventTriggerPrimTrackTrueMultHist1,minMCEventTriggerPrimTrackTrueMultHist1,maxMCEventTriggerPrimTrackTrueMultHist1);
471   fMCTriggerPrimTrackTrueMultHist1->SetBinEdges(0,binsPt);
472   fMCTriggerPrimTrackTrueMultHist1->SetBinEdges(1,binsEta);
473   fMCTriggerPrimTrackTrueMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
474   fMCTriggerPrimTrackTrueMultHist1->GetAxis(1)->SetTitle("mcEta");
475   fMCTriggerPrimTrackTrueMultHist1->GetAxis(2)->SetTitle("true_mult");
476   fMCTriggerPrimTrackTrueMultHist1->Sumw2();
477
478   Int_t binsMCEventPrimTrackTrueMultHist1[3]={ptNbins,etaNbins,150};
479   Double_t minMCEventPrimTrackTrueMultHist1[3]={0.,-1.,-0.5}; 
480   Double_t maxMCEventPrimTrackTrueMultHist1[3]={20.,1.,149.5}; 
481   snprintf(name,256,"fMCEventPrimTrackTrueMultHist1");
482   snprintf(title,256,"mcPt:mcEta:true_mult");
483   
484   fMCEventPrimTrackTrueMultHist1 = new THnSparseF(name,title,3,binsMCEventPrimTrackTrueMultHist1,minMCEventPrimTrackTrueMultHist1,maxMCEventPrimTrackTrueMultHist1);
485   fMCEventPrimTrackTrueMultHist1->SetBinEdges(0,binsPt);
486   fMCEventPrimTrackTrueMultHist1->SetBinEdges(1,binsEta);
487   fMCEventPrimTrackTrueMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");
488   fMCEventPrimTrackTrueMultHist1->GetAxis(1)->SetTitle("mcEta");
489   fMCEventPrimTrackTrueMultHist1->GetAxis(2)->SetTitle("true_mult");
490   fMCEventPrimTrackTrueMultHist1->Sumw2();
491
492   //
493   // mcPT vs multiplicity vs true multiplicity
494   //
495
496   Int_t binsMCAllPrimTrackTrueMultHist2[3]={ptNbins,150,150};
497   Double_t minMCAllPrimTrackTrueMultHist2[3]={0.,-0.5,-0.5}; 
498   Double_t maxMCAllPrimTrackTrueMultHist2[3]={20.,149.5,149.5}; 
499   snprintf(name,256,"fMCAllPrimTrackTrueMultHist2");
500   snprintf(title,256,"mcPt:mult:true_mult");
501   
502   fMCAllPrimTrackTrueMultHist2 = new THnSparseF(name,title,3,binsMCAllPrimTrackTrueMultHist2,minMCAllPrimTrackTrueMultHist2,maxMCAllPrimTrackTrueMultHist2);
503   fMCAllPrimTrackTrueMultHist2->SetBinEdges(0,binsPt);
504   fMCAllPrimTrackTrueMultHist2->GetAxis(0)->SetTitle("mcPt (GeV/c)");
505   fMCAllPrimTrackTrueMultHist2->GetAxis(1)->SetTitle("mult");
506   fMCAllPrimTrackTrueMultHist2->GetAxis(2)->SetTitle("true_mult");
507   fMCAllPrimTrackTrueMultHist2->Sumw2();
508
509   Int_t binsMCNDEventAllPrimTrackTrueMultHist2[3]={ptNbins,150,150};
510   Double_t minMCNDEventAllPrimTrackTrueMultHist2[3]={0.,-0.5,-0.5}; 
511   Double_t maxMCNDEventAllPrimTrackTrueMultHist2[3]={20.,149.5,149.5}; 
512   snprintf(name,256,"fMCNDEventAllPrimTrackTrueMultHist2");
513   snprintf(title,256,"mcPt:mult:true_mult");
514   
515   fMCNDEventAllPrimTrackTrueMultHist2 = new THnSparseF(name,title,3,binsMCNDEventAllPrimTrackTrueMultHist2,minMCNDEventAllPrimTrackTrueMultHist2,maxMCNDEventAllPrimTrackTrueMultHist2);
516   fMCNDEventAllPrimTrackTrueMultHist2->SetBinEdges(0,binsPt);
517   fMCNDEventAllPrimTrackTrueMultHist2->GetAxis(0)->SetTitle("mcPt (GeV/c)");
518   fMCNDEventAllPrimTrackTrueMultHist2->GetAxis(1)->SetTitle("mult");
519   fMCNDEventAllPrimTrackTrueMultHist2->GetAxis(2)->SetTitle("true_mult");
520   fMCNDEventAllPrimTrackTrueMultHist2->Sumw2();
521
522   Int_t binsMCNSDEventAllPrimTrackTrueMultHist2[3]={ptNbins,150,150};
523   Double_t minMCNSDEventAllPrimTrackTrueMultHist2[3]={0.,-0.5,-0.5}; 
524   Double_t maxMCNSDEventAllPrimTrackTrueMultHist2[3]={20.,149.5,149.5}; 
525   snprintf(name,256,"fMCNSDEventAllPrimTrackTrueMultHist2");
526   snprintf(title,256,"mcPt:mult:true_mult");
527   
528   fMCNSDEventAllPrimTrackTrueMultHist2 = new THnSparseF(name,title,3,binsMCNSDEventAllPrimTrackTrueMultHist2,minMCNSDEventAllPrimTrackTrueMultHist2,maxMCNSDEventAllPrimTrackTrueMultHist2);
529   fMCNSDEventAllPrimTrackTrueMultHist2->SetBinEdges(0,binsPt);
530   fMCNSDEventAllPrimTrackTrueMultHist2->GetAxis(0)->SetTitle("mcPt (GeV/c)");
531   fMCNSDEventAllPrimTrackTrueMultHist2->GetAxis(1)->SetTitle("mult");
532   fMCNSDEventAllPrimTrackTrueMultHist2->GetAxis(2)->SetTitle("true_mult");
533   fMCNSDEventAllPrimTrackTrueMultHist2->Sumw2();
534
535   Int_t binsMCEventTriggerPrimTrackTrueMultHist2[3]={ptNbins,150,150};
536   Double_t minMCEventTriggerPrimTrackTrueMultHist2[3]={0.,-0.5,-0.5}; 
537   Double_t maxMCEventTriggerPrimTrackTrueMultHist2[3]={20.,149.5,149.5}; 
538   snprintf(name,256,"fMCTriggerPrimTrackTrueMultHist2");
539   snprintf(title,256,"mcPt:mult:true_mult");
540   
541   fMCTriggerPrimTrackTrueMultHist2 = new THnSparseF(name,title,3,binsMCEventTriggerPrimTrackTrueMultHist2,minMCEventTriggerPrimTrackTrueMultHist2,maxMCEventTriggerPrimTrackTrueMultHist2);
542   fMCTriggerPrimTrackTrueMultHist2->SetBinEdges(0,binsPt);
543   fMCTriggerPrimTrackTrueMultHist2->GetAxis(0)->SetTitle("mcPt (GeV/c)");
544   fMCTriggerPrimTrackTrueMultHist2->GetAxis(1)->SetTitle("mult");
545   fMCTriggerPrimTrackTrueMultHist2->GetAxis(2)->SetTitle("true_mult");
546   fMCTriggerPrimTrackTrueMultHist2->Sumw2();
547
548   Int_t binsMCEventPrimTrackTrueMultHist2[3]={ptNbins,150,150};
549   Double_t minMCEventPrimTrackTrueMultHist2[3]={0.,-0.5,-0.5}; 
550   Double_t maxMCEventPrimTrackTrueMultHist2[3]={20.,149.5,149.5}; 
551   snprintf(name,256,"fMCEventPrimTrackTrueMultHist2");
552   snprintf(title,256,"mcPt:mult:true_mult");
553   
554   fMCEventPrimTrackTrueMultHist2 = new THnSparseF(name,title,3,binsMCEventPrimTrackTrueMultHist2,minMCEventPrimTrackTrueMultHist2,maxMCEventPrimTrackTrueMultHist2);
555   fMCEventPrimTrackTrueMultHist2->SetBinEdges(0,binsPt);
556   fMCEventPrimTrackTrueMultHist2->GetAxis(0)->SetTitle("mcPt (GeV/c)");
557   fMCEventPrimTrackTrueMultHist2->GetAxis(1)->SetTitle("mult");
558   fMCEventPrimTrackTrueMultHist2->GetAxis(2)->SetTitle("true_mult");
559   fMCEventPrimTrackTrueMultHist2->Sumw2();
560
561
562   //
563   // mean pt
564   //
565   Int_t binsMCAllPrimTrackMeanPtTrueMult1[2]={100,150};
566   Double_t minMCAllPrimTrackMeanPtTrueMult1[2]={0.,-0.5}; 
567   Double_t maxMCAllPrimTrackMeanPtTrueMult1[2]={10.,149.5}; 
568   snprintf(name,256,"fMCAllPrimTrackMeanPtTrueMult1");
569   snprintf(title,256,"event <mcPt>:true_mult");
570   
571   fMCAllPrimTrackMeanPtTrueMult1 = new THnSparseF(name,title,2,binsMCAllPrimTrackMeanPtTrueMult1,minMCAllPrimTrackMeanPtTrueMult1,maxMCAllPrimTrackMeanPtTrueMult1);
572   fMCAllPrimTrackMeanPtTrueMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
573   fMCAllPrimTrackMeanPtTrueMult1->GetAxis(1)->SetTitle("true_mult");
574   fMCAllPrimTrackMeanPtTrueMult1->Sumw2();
575
576   Int_t binsMCAllPrimTrackMeanPtMult1[2]={100,150};
577   Double_t minMCAllPrimTrackMeanPtMult1[2]={0.,-0.5}; 
578   Double_t maxMCAllPrimTrackMeanPtMult1[2]={10.,149.5}; 
579   snprintf(name,256,"fMCAllPrimTrackMeanPtMult1");
580   snprintf(title,256,"event <mcPt>:mult");
581   
582   fMCAllPrimTrackMeanPtMult1 = new THnSparseF(name,title,2,binsMCAllPrimTrackMeanPtMult1,minMCAllPrimTrackMeanPtMult1,maxMCAllPrimTrackMeanPtMult1);
583   fMCAllPrimTrackMeanPtMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
584   fMCAllPrimTrackMeanPtMult1->GetAxis(1)->SetTitle("multiplicity");
585   fMCAllPrimTrackMeanPtMult1->Sumw2();
586
587   //
588   Int_t binsMCNDEventAllPrimTrackMeanPtTrueMult1[2]={100,150};
589   Double_t minMCNDEventAllPrimTrackMeanPtTrueMult1[2]={0.,-0.5}; 
590   Double_t maxMCNDEventAllPrimTrackMeanPtTrueMult1[2]={10.,149.5}; 
591   snprintf(name,256,"fMCNDEventAllPrimTrackMeanPtTrueMult1");
592   snprintf(title,256,"event <mcPt>:true_mult");
593   
594   fMCNDEventAllPrimTrackMeanPtTrueMult1 = new THnSparseF(name,title,2,binsMCNDEventAllPrimTrackMeanPtTrueMult1,minMCNDEventAllPrimTrackMeanPtTrueMult1,maxMCNDEventAllPrimTrackMeanPtTrueMult1);
595   fMCNDEventAllPrimTrackMeanPtTrueMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
596   fMCNDEventAllPrimTrackMeanPtTrueMult1->GetAxis(1)->SetTitle("true_mult");
597   fMCNDEventAllPrimTrackMeanPtTrueMult1->Sumw2();
598
599   Int_t binsMCNDEventAllPrimTrackMeanPtMult1[2]={100,150};
600   Double_t minMCNDEventAllPrimTrackMeanPtMult1[2]={0.,-0.5}; 
601   Double_t maxMCNDEventAllPrimTrackMeanPtMult1[2]={10.,149.5}; 
602   snprintf(name,256,"fMCNDEventAllPrimTrackMeanPtMult1");
603   snprintf(title,256,"event <mcPt>:mult");
604   
605   fMCNDEventAllPrimTrackMeanPtMult1 = new THnSparseF(name,title,2,binsMCNDEventAllPrimTrackMeanPtMult1,minMCNDEventAllPrimTrackMeanPtMult1,maxMCNDEventAllPrimTrackMeanPtMult1);
606   fMCNDEventAllPrimTrackMeanPtMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
607   fMCNDEventAllPrimTrackMeanPtMult1->GetAxis(1)->SetTitle("multiplicity");
608   fMCNDEventAllPrimTrackMeanPtMult1->Sumw2();
609
610   //
611   Int_t binsMCNSDEventAllPrimTrackMeanPtTrueMult1[2]={100,150};
612   Double_t minMCNSDEventAllPrimTrackMeanPtTrueMult1[2]={0.,-0.5}; 
613   Double_t maxMCNSDEventAllPrimTrackMeanPtTrueMult1[2]={10.,149.5}; 
614   snprintf(name,256,"fMCNSDEventAllPrimTrackMeanPtTrueMult1");
615   snprintf(title,256,"event <mcPt>:true_mult");
616   
617   fMCNSDEventAllPrimTrackMeanPtTrueMult1 = new THnSparseF(name,title,2,binsMCNSDEventAllPrimTrackMeanPtTrueMult1,minMCNSDEventAllPrimTrackMeanPtTrueMult1,maxMCNSDEventAllPrimTrackMeanPtTrueMult1);
618   fMCNSDEventAllPrimTrackMeanPtTrueMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
619   fMCNSDEventAllPrimTrackMeanPtTrueMult1->GetAxis(1)->SetTitle("true_mult");
620   fMCNSDEventAllPrimTrackMeanPtTrueMult1->Sumw2();
621
622   Int_t binsMCNSDEventAllPrimTrackMeanPtMult1[2]={100,150};
623   Double_t minMCNSDEventAllPrimTrackMeanPtMult1[2]={0.,-0.5}; 
624   Double_t maxMCNSDEventAllPrimTrackMeanPtMult1[2]={10.,149.5}; 
625   snprintf(name,256,"fMCNSDEventAllPrimTrackMeanPtMult1");
626   snprintf(title,256,"event <mcPt>:mult");
627   
628   fMCNSDEventAllPrimTrackMeanPtMult1 = new THnSparseF(name,title,2,binsMCNSDEventAllPrimTrackMeanPtMult1,minMCNSDEventAllPrimTrackMeanPtMult1,maxMCNSDEventAllPrimTrackMeanPtMult1);
629   fMCNSDEventAllPrimTrackMeanPtMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
630   fMCNSDEventAllPrimTrackMeanPtMult1->GetAxis(1)->SetTitle("multiplicity");
631   fMCNSDEventAllPrimTrackMeanPtMult1->Sumw2();
632
633   //
634   Int_t binsMCTriggerPrimTrackMeanPtTrueMult1[2]={100,150};
635   Double_t minMCTriggerPrimTrackMeanPtTrueMult1[2]={0.,-0.5}; 
636   Double_t maxMCTriggerPrimTrackMeanPtTrueMult1[2]={10.,149.5}; 
637   snprintf(name,256,"fMCTriggerPrimTrackMeanPtTrueMult1");
638   snprintf(title,256,"event <mcPt>:true_mult");
639   
640   fMCTriggerPrimTrackMeanPtTrueMult1 = new THnSparseF(name,title,2,binsMCTriggerPrimTrackMeanPtTrueMult1,minMCTriggerPrimTrackMeanPtTrueMult1,maxMCTriggerPrimTrackMeanPtTrueMult1);
641   fMCTriggerPrimTrackMeanPtTrueMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
642   fMCTriggerPrimTrackMeanPtTrueMult1->GetAxis(1)->SetTitle("true_mult");
643   fMCTriggerPrimTrackMeanPtTrueMult1->Sumw2();
644
645   Int_t binsMCTriggerPrimTrackMeanPtMult1[2]={100,150};
646   Double_t minMCTriggerPrimTrackMeanPtMult1[2]={0.,-0.5}; 
647   Double_t maxMCTriggerPrimTrackMeanPtMult1[2]={10.,149.5}; 
648   snprintf(name,256,"fMCTriggerPrimTrackMeanPtMult1");
649   snprintf(title,256,"event <mcPt>:mult");
650   
651   fMCTriggerPrimTrackMeanPtMult1 = new THnSparseF(name,title,2,binsMCTriggerPrimTrackMeanPtMult1,minMCTriggerPrimTrackMeanPtMult1,maxMCTriggerPrimTrackMeanPtMult1);
652   fMCTriggerPrimTrackMeanPtMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
653   fMCTriggerPrimTrackMeanPtMult1->GetAxis(1)->SetTitle("multiplicity");
654   fMCTriggerPrimTrackMeanPtMult1->Sumw2();
655
656   //
657   Int_t binsMCEventPrimTrackMeanPtTrueMult1[2]={100,150};
658   Double_t minMCEventPrimTrackMeanPtTrueMult1[2]={0.,-0.5}; 
659   Double_t maxMCEventPrimTrackMeanPtTrueMult1[2]={10.,149.5}; 
660   snprintf(name,256,"fMCEventPrimTrackMeanPtTrueMult1");
661   snprintf(title,256,"event <mcPt>:true_mult");
662   
663   fMCEventPrimTrackMeanPtTrueMult1 = new THnSparseF(name,title,2,binsMCEventPrimTrackMeanPtTrueMult1,minMCEventPrimTrackMeanPtTrueMult1,maxMCEventPrimTrackMeanPtTrueMult1);
664   fMCEventPrimTrackMeanPtTrueMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
665   fMCEventPrimTrackMeanPtTrueMult1->GetAxis(1)->SetTitle("true_mult");
666   fMCEventPrimTrackMeanPtTrueMult1->Sumw2();
667
668   Int_t binsMCEventPrimTrackMeanPtMult1[2]={100,150};
669   Double_t minMCEventPrimTrackMeanPtMult1[2]={0.,-0.5}; 
670   Double_t maxMCEventPrimTrackMeanPtMult1[2]={10.,149.5}; 
671   snprintf(name,256,"fMCEventPrimTrackMeanPtMult1");
672   snprintf(title,256,"event <mcPt>:mult");
673   
674   fMCEventPrimTrackMeanPtMult1 = new THnSparseF(name,title,2,binsMCEventPrimTrackMeanPtMult1,minMCEventPrimTrackMeanPtMult1,maxMCEventPrimTrackMeanPtMult1);
675   fMCEventPrimTrackMeanPtMult1->GetAxis(0)->SetTitle("<mcPt> (GeV/c)");
676   fMCEventPrimTrackMeanPtMult1->GetAxis(1)->SetTitle("multiplicity");
677   fMCEventPrimTrackMeanPtMult1->Sumw2();
678
679
680
681
682
683
684   //
685   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) 
686   {
687     // THnSparse track histograms
688     //
689     Int_t binsRecTrackHist1[3]={ptNbins,etaNbins,90};
690     Double_t minRecTrackHist1[3]={0.,-1.,0.}; 
691     Double_t maxRecTrackHist1[3]={10.,1.,2.*TMath::Pi()};
692     snprintf(name,256,"fRecTrackHist1_%d",i);
693     snprintf(title,256,"Pt:Eta:Phi");
694   
695     fRecTrackHist1[i] = new THnSparseF(name,title,3,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1);
696     fRecTrackHist1[i]->SetBinEdges(0,binsPt);
697     fRecTrackHist1[i]->SetBinEdges(1,binsEta);
698     fRecTrackHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
699     fRecTrackHist1[i]->GetAxis(1)->SetTitle("Eta");
700     fRecTrackHist1[i]->GetAxis(2)->SetTitle("Phi (rad)");
701     fRecTrackHist1[i]->Sumw2();
702   }
703
704   //
705   Int_t binsCorrRecTrackMultHist1[3]={ptNbins,etaNbins,150};
706   Double_t minCorrRecTrackMultHist1[3]={0.,-1.,-0.5}; 
707   Double_t maxCorrRecTrackMultHist1[3]={20.,1.,149.5};
708
709   Int_t binsCorrRecTrackTrueMultHist1[3]={ptNbins,etaNbins,150};
710   Double_t minCorrRecTrackTrueMultHist1[3]={0.,-1.,-0.5}; 
711   Double_t maxCorrRecTrackTrueMultHist1[3]={20.,1.,149.5};
712
713   Int_t binsCorrRecTrackTrueMultHist2[3]={ptNbins,150,150};
714   Double_t minCorrRecTrackTrueMultHist2[3]={0.,-0.5,-0.5}; 
715   Double_t maxCorrRecTrackTrueMultHist2[3]={20.,149.5,149.5};
716
717   //
718   Int_t binsCorrRecTrackMeanPtMultHist1[2]={100,150};
719   Double_t minCorrRecTrackMeanPtMultHist1[2]={0.,-0.5}; 
720   Double_t maxCorrRecTrackMeanPtMultHist1[2]={10.,149.5};
721
722   Int_t binsCorrRecTrackMeanPtTrueMultHist1[2]={100,150};
723   Double_t minCorrRecTrackMeanPtTrueMultHist1[2]={0.,-0.5}; 
724   Double_t maxCorrRecTrackMeanPtTrueMultHist1[2]={10.,149.5};
725
726   Int_t binsCorrRecTrackPt1[1]={200};
727   Double_t minCorrRecTrackPt1[1]={0.}; 
728   Double_t maxCorrRecTrackPt1[1]={10.};
729
730   for(Int_t i=0; i<8; i++) 
731   {
732     // THnSparse track histograms
733     snprintf(name,256,"fCorrRecTrackMultHist1_%d",i);
734     snprintf(title,256,"Pt:Eta:mult");
735     fCorrRecTrackMultHist1[i] = new THnSparseF(name,title,3,binsCorrRecTrackMultHist1,minCorrRecTrackMultHist1,maxCorrRecTrackMultHist1);
736     fCorrRecTrackMultHist1[i]->SetBinEdges(0,binsPt);
737     fCorrRecTrackMultHist1[i]->SetBinEdges(1,binsEta);
738     fCorrRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
739     fCorrRecTrackMultHist1[i]->GetAxis(1)->SetTitle("Eta");
740     fCorrRecTrackMultHist1[i]->GetAxis(2)->SetTitle("multiplicity");
741     fCorrRecTrackMultHist1[i]->Sumw2();
742
743     // THnSparse track histograms
744     snprintf(name,256,"fCorrRecTrackTrueMultHist1_%d",i);
745     snprintf(title,256,"Pt:Eta:true_mult");
746     fCorrRecTrackTrueMultHist1[i] = new THnSparseF(name,title,3,binsCorrRecTrackTrueMultHist1,minCorrRecTrackTrueMultHist1,maxCorrRecTrackTrueMultHist1);
747     fCorrRecTrackTrueMultHist1[i]->SetBinEdges(0,binsPt);
748     fCorrRecTrackTrueMultHist1[i]->SetBinEdges(1,binsEta);
749     fCorrRecTrackTrueMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
750     fCorrRecTrackTrueMultHist1[i]->GetAxis(1)->SetTitle("Eta");
751     fCorrRecTrackTrueMultHist1[i]->GetAxis(2)->SetTitle("true multiplicity");
752     fCorrRecTrackTrueMultHist1[i]->Sumw2();
753
754     //
755     snprintf(name,256,"fCorrRecTrackTrueMultHist2_%d",i);
756     snprintf(title,256,"Pt:mult:true_mult");
757     fCorrRecTrackTrueMultHist2[i] = new THnSparseF(name,title,3,binsCorrRecTrackTrueMultHist2,minCorrRecTrackTrueMultHist2,maxCorrRecTrackTrueMultHist2);
758     fCorrRecTrackTrueMultHist2[i]->SetBinEdges(0,binsPt);
759     fCorrRecTrackTrueMultHist2[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
760     fCorrRecTrackTrueMultHist2[i]->GetAxis(1)->SetTitle("multiplicity");
761     fCorrRecTrackTrueMultHist2[i]->GetAxis(2)->SetTitle("true multiplicity");
762     fCorrRecTrackTrueMultHist2[i]->Sumw2();
763
764     // THnSparse track histograms
765     snprintf(name,256,"fCorrRecTrackMeanPtMultHist1_%d",i);
766     snprintf(title,256,"<Pt>:mult");
767     fCorrRecTrackMeanPtMultHist1[i] = new THnSparseF(name,title,2,binsCorrRecTrackMeanPtMultHist1,minCorrRecTrackMeanPtMultHist1,maxCorrRecTrackMeanPtMultHist1);
768     fCorrRecTrackMeanPtMultHist1[i]->GetAxis(0)->SetTitle("<Pt> (GeV/c)");
769     fCorrRecTrackMeanPtMultHist1[i]->GetAxis(1)->SetTitle("multiplicity");
770     fCorrRecTrackMeanPtMultHist1[i]->Sumw2();
771
772     // THnSparse track histograms
773     snprintf(name,256,"fCorrRecTrackMeanPtTrueMultHist1_%d",i);
774     snprintf(title,256,"<Pt>:true_mult");
775     fCorrRecTrackMeanPtTrueMultHist1[i] = new THnSparseF(name,title,2,binsCorrRecTrackMeanPtTrueMultHist1,minCorrRecTrackMeanPtTrueMultHist1,maxCorrRecTrackMeanPtTrueMultHist1);
776     fCorrRecTrackMeanPtTrueMultHist1[i]->GetAxis(0)->SetTitle("<Pt> (GeV/c)");
777     fCorrRecTrackMeanPtTrueMultHist1[i]->GetAxis(1)->SetTitle("true multiplicity");
778     fCorrRecTrackMeanPtTrueMultHist1[i]->Sumw2();
779
780     snprintf(name,256,"fCorrRecTrackPt1_%d",i);
781     snprintf(title,256,"pt small bining");
782     fCorrRecTrackPt1[i] = new THnSparseF(name,title,1,binsCorrRecTrackPt1,minCorrRecTrackPt1,maxCorrRecTrackPt1);
783     fCorrRecTrackPt1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");
784     fCorrRecTrackPt1[i]->Sumw2();
785
786   }
787
788   Int_t binsEventMatrix[2]={zvNbins,150};
789   Double_t minEventMatrix[2]={-25.,-0.5};
790   Double_t maxEventMatrix[2]={25.,149.5};
791
792   fMCAllEventMultHist1 = new THnSparseF("fMCAllEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
793   fMCAllEventMultHist1->SetBinEdges(0,binsZv);
794   fMCAllEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");
795   fMCAllEventMultHist1->GetAxis(1)->SetTitle("multiplicity");
796   fMCAllEventMultHist1->Sumw2();
797
798   fMCAllNDEventMultHist1 = new THnSparseF("fMCAllNDEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
799   fMCAllNDEventMultHist1->SetBinEdges(0,binsZv);
800   fMCAllNDEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");
801   fMCAllNDEventMultHist1->GetAxis(1)->SetTitle("multiplicity");
802   fMCAllNDEventMultHist1->Sumw2();
803
804   fMCAllNSDEventMultHist1 = new THnSparseF("fMCAllNSDEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
805   fMCAllNSDEventMultHist1->SetBinEdges(0,binsZv);
806   fMCAllNSDEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");
807   fMCAllNSDEventMultHist1->GetAxis(1)->SetTitle("multiplicity");
808   fMCAllNSDEventMultHist1->Sumw2();
809
810   fMCTriggerMultHist1 = new THnSparseF("fMCTriggerMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
811   fMCTriggerMultHist1->SetBinEdges(0,binsZv);
812   fMCTriggerMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");
813   fMCTriggerMultHist1->GetAxis(1)->SetTitle("multiplicity");
814   fMCTriggerMultHist1->Sumw2();
815
816   fMCEventMultHist1 = new THnSparseF("fMCEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
817   fMCEventMultHist1->SetBinEdges(0,binsZv);
818   fMCEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");
819   fMCEventMultHist1->GetAxis(1)->SetTitle("multiplicity");
820   fMCEventMultHist1->Sumw2();
821
822   for(Int_t i=0; i<5; i++) 
823   {
824     // event corrected histograms
825     snprintf(name,256,"fCorrRecEventHist1_%d",i);
826     snprintf(title,256,"mcZv:mult");
827     fCorrRecEventHist1[i] = new THnSparseF("fCorrRecEventHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
828     fCorrRecEventHist1[i]->SetBinEdges(0,binsZv);
829     fCorrRecEventHist1[i]->GetAxis(0)->SetTitle("Zv (cm)");
830     fCorrRecEventHist1[i]->GetAxis(1)->SetTitle("multiplicity");
831     fCorrRecEventHist1[i]->Sumw2();
832
833     // empty event corrected histograms
834     snprintf(name,256,"fCorrRecEventHist2_%d",i);
835     snprintf(title,256,"mcZv:mult");
836     fCorrRecEventHist2[i] = new THnSparseF("fCorrRecEventHist2","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);
837     fCorrRecEventHist2[i]->SetBinEdges(0,binsZv);
838     fCorrRecEventHist2[i]->GetAxis(0)->SetTitle("Zv (cm)");
839     fCorrRecEventHist2[i]->GetAxis(1)->SetTitle("multiplicity");
840     fCorrRecEventHist2[i]->Sumw2();
841   }
842   
843   //
844   // cosmics histo
845   //
846   Int_t binsCosmicsHisto[3]=  {151, 300, ptNbins};
847   Double_t minCosmicsHisto[3]={-1.5, -2.*TMath::Pi(), 0.0}; 
848   Double_t maxCosmicsHisto[3]={ 1.5, 2.*TMath::Pi(), 16.0}; 
849   snprintf(name,256,"fCosmicsHisto");
850   snprintf(title,256,"deta:dphi:pt");
851   
852   fCosmicsHisto = new THnSparseF(name,title,3,binsCosmicsHisto,minCosmicsHisto,maxCosmicsHisto);
853   fCosmicsHisto->SetBinEdges(2,binsPt);
854   fCosmicsHisto->GetAxis(0)->SetTitle("#Delta#eta");
855   fCosmicsHisto->GetAxis(1)->SetTitle("#Delta#phi (rad)");
856   fCosmicsHisto->GetAxis(2)->SetTitle("pt (GV/c)");
857   fCosmicsHisto->Sumw2();
858
859   //
860   Int_t binsEventCount[3]={2,2,2};
861   Double_t minEventCount[3]={0,0,0}; 
862   Double_t maxEventCount[3]={2,2,2}; 
863   fEventCount = new THnSparseF("fEventCount","trig vs trig+vertex",3,binsEventCount,minEventCount,maxEventCount);
864   fEventCount->GetAxis(0)->SetTitle("trig");
865   fEventCount->GetAxis(1)->SetTitle("trig+vert");
866   fEventCount->GetAxis(2)->SetTitle("selected");
867   fEventCount->Sumw2();
868
869
870   // init output folder
871   fCorrectionFolder = CreateFolder("folderdNdPt","Correction dNdPt Folder");
872
873   // init correction matrices
874   TFile *file = (TFile *)TFile::Open(fCorrMatrixFileName.Data()); 
875   if(!file) { 
876     AliDebug(AliLog::kError, "file with efficiency matrices not available");
877     printf("file with efficiency matrices not available \n");
878   } else {
879     TFolder *folder = (TFolder *)file->FindObjectAny("folderdNdPt");
880     if(!folder) { 
881       AliDebug(AliLog::kError, "file without folderdNdPt");
882       printf("file without folderdNdPt \n");
883     } else {
884       // rec. event mult vs true multiplicity 
885       fEventMultCorrelationMatrix = (THnSparseF*)folder->FindObject("event_mult_correlation_matrix");
886       if(!fEventMultCorrelationMatrix) {
887          Printf("No %s matrix \n", "event_mult_correlation_matrix");
888          return;
889       }
890
891       //
892       // event level corrections (zv,mult_MB)
893       //
894  
895       // trigger bias correction (MBtoND) 
896       fCorrTriggerMBtoNDEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_trig_MBtoND_corr_matrix");
897       if(!fCorrTriggerMBtoNDEventMatrix) {
898          Printf("No %s matrix \n", "zv_mult_trig_MBtoND_corr_matrix");
899          return;
900       }
901
902       // trigger bias correction (MBtoNSD)
903       fCorrTriggerMBtoNSDEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_trig_MBtoNSD_corr_matrix");
904       if(!fCorrTriggerMBtoNSDEventMatrix) {
905          Printf("No %s matrix \n", "zv_mult_trig_MBtoNSD_corr_matrix");
906          return;
907       }
908
909       // trigger bias correction (MBtoInel)
910       fCorrTriggerMBtoInelEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_trig_MBtoInel_corr_matrix");
911       if(!fCorrTriggerMBtoInelEventMatrix) {
912          Printf("No %s matrix \n", "zv_mult_trig_MBtoInel_corr_matrix"); 
913          return;
914       }
915      
916       // vertex reconstruction efficiency correction
917       fCorrEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_event_corr_matrix");
918       if(!fCorrEventMatrix) {
919          Printf("No %s matrix \n", "zv_mult_event_corr_matrix");
920          return;
921       }
922
923       //
924       // histogram needed for empty events corrections
925       //
926       fZvNorm = (TH1D*)folder->FindObject("zv_distribution_norm");
927       if(!fZvNorm) {
928          Printf("No %s matrix \n", "fZvNorm");
929          return;
930       }
931
932       fZvEmptyEventsNorm = (TH1D*)folder->FindObject("zv_empty_events_norm");
933       if(!fZvEmptyEventsNorm) {
934          Printf("No %s matrix \n", "fZvEmptyEventsNorm");
935          return;
936       }
937
938       fLHCBin0Background = (TH1D*)folder->FindObject("hLHCBin0Background");
939       if(!fLHCBin0Background) {
940          Printf("No %s matrix \n", "fLHCBin0Background");
941          return;
942       }
943
944       //
945       // track-event level corrections (zv,pt,eta)
946       //
947
948       // trigger bias correction (MBtoND) 
949       fCorrTriggerMBtoNDTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_trig_MBtoND_corr_matrix");
950       if(!fCorrTriggerMBtoNDTrackEventMatrix) {
951          Printf("No %s matrix \n", "zv_pt_eta_track_trig_MBtoND_corr_matrix");
952          return;
953       }
954
955       // trigger bias correction (MBtoNSD)
956       fCorrTriggerMBtoNSDTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_trig_MBtoNSD_corr_matrix");
957       if(!fCorrTriggerMBtoNSDTrackEventMatrix) {
958          Printf("No %s matrix \n", "zv_pt_eta_track_trig_MBtoNSD_corr_matrix");
959          return;
960       }
961
962       // trigger bias correction (MBtoInel) 
963       fCorrTriggerMBtoInelTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_trig_MBtoInel_corr_matrix");
964       if(!fCorrTriggerMBtoInelTrackEventMatrix) {
965          Printf("No %s matrix \n", "zv_pt_eta_track_trig_MBtoInel_corr_matrix");
966          return;
967       }
968     
969       // vertex reconstruction efficiency correction (zv,pt,eta)
970       fCorrTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_event_corr_matrix");
971       if(!fCorrTrackEventMatrix) {
972          Printf("No %s matrix \n", "zv_pt_eta_track_event_corr_matrix");
973          return;
974       }
975
976       // track reconstruction efficiency correction (zv,pt,eta)
977       fCorrTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_corr_matrix");
978       if(!fCorrTrackMatrix) {
979          Printf("No %s matrix \n", "zv_pt_eta_track_corr_matrix");
980          return;
981       }
982
983       // high pt track reconstruction efficiency correction (zv,pt,eta)
984       fCorrHighPtTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_highPt_track_corr_matrix");
985       if(!fCorrHighPtTrackMatrix) {
986          Printf("No %s matrix \n", "zv_pt_eta_highPt_track_corr_matrix");
987          return;
988       }
989
990       // secondary tracks contamination correction (zv,pt,eta)
991       fContTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_cont_matrix");
992       if(!fContTrackMatrix) {
993          Printf("No %s matrix \n", "zv_pt_eta_track_cont_matrix");
994          return;
995       }
996
997       // multiply reconstructed tracks correction
998       fContMultTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_mult_track_cont_matrix");
999       if(!fContMultTrackMatrix) {
1000          Printf("No %s matrix \n", "zv_pt_eta_mult_track_cont_matrix");
1001          return;
1002       }
1003     }
1004   }
1005
1006 }
1007
1008 //_____________________________________________________________________________
1009 void AlidNdPtCorrection::Process(AliESDEvent *esdEvent, AliMCEvent *mcEvent)
1010 {
1011   //
1012   // Process real and/or simulated events
1013   //
1014   if(!esdEvent) {
1015     AliDebug(AliLog::kError, "esdEvent not available");
1016     return;
1017   }
1018
1019   // get selection cuts
1020   AlidNdPtEventCuts *evtCuts = GetEventCuts(); 
1021   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); 
1022   AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); 
1023
1024   if(!evtCuts || !accCuts  || !esdTrackCuts) {
1025     AliDebug(AliLog::kError, "cuts not available");
1026     return;
1027   }
1028
1029
1030   // trigger selection
1031   Bool_t isEventTriggered = kTRUE;
1032   AliPhysicsSelection *trigSel = NULL;
1033   AliTriggerAnalysis *trigAna = NULL; // needed for andV0
1034
1035   if(evtCuts->IsTriggerRequired())  
1036   {
1037     //
1038     trigSel = GetPhysicsTriggerSelection();
1039     if(!trigSel) {
1040       AliDebug(AliLog::kError, "cannot get trigSel");
1041       return;
1042     }
1043     
1044     if(IsUseMCInfo()) 
1045     { 
1046       trigSel->SetAnalyzeMC();
1047
1048       if(GetParticleMode() == AlidNdPtHelper::kVZEROCase1)
1049       {
1050         // check V0 systematics (case1)
1051         // Initialization done in the macro
1052         trigAna = trigSel->GetTriggerAnalysis();
1053         if(!trigAna) 
1054           return;
1055
1056         //trigAna->SetV0HwPars(15, 61.5, 86.5);
1057         //trigAna->SetV0AdcThr(15);
1058
1059         isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);
1060         //printf("MB1 & kVZEROCase1  %d \n",isEventTriggered);
1061         //isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1062         
1063         if(GetTrigger() == AliTriggerAnalysis::kV0AND) 
1064         {
1065           isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1066           //printf("V0AND %d \n",isEventTriggered);
1067         }
1068       }
1069       else if(GetParticleMode() == AlidNdPtHelper::kVZEROCase2)
1070       {
1071         // check V0 systematics (case2 only in MC)
1072         // Initialization done in the macro
1073
1074         trigAna = trigSel->GetTriggerAnalysis();
1075         if(!trigAna) 
1076           return;
1077
1078         //trigAna->SetV0HwPars(0, 0, 125);
1079         //trigAna->SetV0AdcThr(0);
1080
1081         isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);
1082         //isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1083         
1084         if(GetTrigger() == AliTriggerAnalysis::kV0AND) 
1085         {
1086           isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1087           //printf("V0AND %d \n",isEventTriggered);
1088         }
1089       }
1090       else {
1091         isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);
1092         //printf("MB1 %d \n",isEventTriggered);
1093         
1094         if(GetTrigger() == AliTriggerAnalysis::kV0AND) 
1095         {
1096           trigAna = trigSel->GetTriggerAnalysis();
1097           if(!trigAna) 
1098             return;
1099
1100           isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1101           //printf("V0AND %d \n",isEventTriggered);
1102         }
1103       }
1104     }
1105     else {
1106       //
1107       // 0-multiplicity bin for LHC background correction
1108       //
1109       if( GetAnalysisMode() == AlidNdPtHelper::kTPCITS || 
1110           GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtx || 
1111           GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtxUpdate || 
1112           GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtx || 
1113           GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt ) 
1114       {
1115         trigSel->SetBin0CallbackViaPointer(&AlidNdPtAnalysis::IsBinZeroTrackSPDvtx);
1116       } else {
1117         trigSel->SetBin0CallbackViaPointer(&AlidNdPtAnalysis::IsBinZeroSPDvtx);
1118       }
1119
1120       if(GetParticleMode() == AlidNdPtHelper::kVZEROCase1)
1121       {
1122         // check V0 systematics (case1)
1123         // Initialization done in the macro
1124         trigAna = trigSel->GetTriggerAnalysis();
1125         if(!trigAna) 
1126           return;
1127
1128         //trigAna->SetV0HwPars(15, 61.5, 86.5);
1129         //trigAna->SetV0AdcThr(15);
1130
1131         isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);
1132         //isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1133         
1134         if(GetTrigger() == AliTriggerAnalysis::kV0AND) 
1135         {
1136           isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1137           //printf("V0AND %d \n",isEventTriggered);
1138         }
1139       }
1140       else {
1141         isEventTriggered = trigSel->IsCollisionCandidate(esdEvent);
1142         //printf("MB1 %d \n",isEventTriggered);
1143         
1144         if(GetTrigger() == AliTriggerAnalysis::kV0AND) 
1145         {
1146           trigAna = trigSel->GetTriggerAnalysis();
1147           if(!trigAna) 
1148             return;
1149
1150           isEventTriggered = trigAna->IsOfflineTriggerFired(esdEvent, GetTrigger());
1151           //printf("V0AND %d \n",isEventTriggered);
1152         }
1153       }
1154     }
1155   }
1156
1157   // use MC information
1158   AliHeader* header = 0;
1159   AliGenEventHeader* genHeader = 0;
1160   AliStack* stack = 0;
1161   TArrayF vtxMC(3);
1162   AliPWG0Helper::MCProcessType evtType = AliPWG0Helper::kInvalidProcess;
1163   Int_t multMCTrueTracks = 0;
1164
1165   if(IsUseMCInfo())
1166   {
1167     if(!mcEvent) {
1168       AliDebug(AliLog::kError, "mcEvent not available");
1169       return;
1170     }
1171
1172     // get MC event header
1173     header = mcEvent->Header();
1174     if (!header) {
1175       AliDebug(AliLog::kError, "Header not available");
1176       return;
1177     }
1178     // MC particle stack
1179     stack = mcEvent->Stack();
1180     if (!stack) {
1181       AliDebug(AliLog::kError, "Stack not available");
1182       return;
1183     }
1184
1185     // get event type (ND=0x1, DD=0x2, SD=0x4)
1186     evtType = AliPWG0Helper::GetEventProcessType(header);
1187     //Printf("evtType %d \n", evtType);
1188     AliDebug(AliLog::kDebug+1, Form("Found process type %d", evtType));
1189
1190     // get MC vertex
1191     genHeader = header->GenEventHeader();
1192     if (!genHeader) {
1193       AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");
1194       return;
1195     }
1196     genHeader->PrimaryVertex(vtxMC);
1197
1198     // Fill MC event histogram
1199     Double_t vMCEventHist1[3]={vtxMC[0],vtxMC[1],vtxMC[2]};
1200     fMCEventHist1->Fill(vMCEventHist1);
1201
1202     // multipliticy of all MC primary tracks
1203     // in Zvtx, eta ranges
1204     multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);
1205
1206   } // end bUseMC
1207
1208   // get reconstructed vertex  
1209   const AliESDVertex* vtxESD = 0; 
1210   Bool_t isRecVertex = kFALSE;
1211   if(evtCuts->IsRecVertexRequired()) 
1212   {
1213     Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();
1214     Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();
1215     vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints); 
1216     if(!vtxESD) return;
1217     isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, esdEvent->GetPrimaryVertexSPD(), GetAnalysisMode(), kFALSE);
1218   }
1219
1220   if( IsUseMCInfo() && !evtCuts->IsRecVertexRequired() ) {
1221     vtxESD = new AliESDVertex(vtxMC[2],10.,genHeader->NProduced(),"smearMC");
1222     if(!vtxESD) return;
1223     isRecVertex = kTRUE;
1224   }
1225
1226   Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex; 
1227   //printf("isEventOK %d, isEventTriggered %d \n",isEventOK, isEventTriggered);
1228   //printf("GetAnalysisMode() %d \n",GetAnalysisMode());
1229
1230   // vertex contributors
1231   Int_t multMBTracks = 0; 
1232   if(GetAnalysisMode() == AlidNdPtHelper::kTPC) 
1233   {  
1234      if(vtxESD->GetStatus() && isRecVertex)
1235        multMBTracks = AlidNdPtHelper::GetTPCMBTrackMult(esdEvent,evtCuts,accCuts,esdTrackCuts);
1236   } 
1237   else if( GetAnalysisMode() == AlidNdPtHelper::kTPCITS ||  GetAnalysisMode() == AlidNdPtHelper::kTPCSPDvtx || 
1238            GetAnalysisMode()==AlidNdPtHelper::kTPCSPDvtxUpdate || GetAnalysisMode()==AlidNdPtHelper::kTPCITSHybrid ) 
1239   {
1240      const AliMultiplicity* mult = esdEvent->GetMultiplicity();
1241      //if(mult && vtxESD->GetStatus() && isRecVertex)
1242      if(mult)
1243        multMBTracks = mult->GetNumberOfTracklets();
1244     
1245   } 
1246   else if( GetAnalysisMode() == AlidNdPtHelper::kTPCITS || 
1247            GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtx || 
1248            GetAnalysisMode() == AlidNdPtHelper::kTPCTrackSPDvtxUpdate || 
1249            GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtx || 
1250            GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt )
1251   {
1252      if(vtxESD && vtxESD->GetStatus() && isRecVertex)
1253        multMBTracks = vtxESD->GetNContributors();
1254
1255   }
1256   else {
1257     AliDebug(AliLog::kError, Form("Found analysis type %d", GetAnalysisMode()));
1258     return; 
1259   }
1260
1261   Bool_t isEventSelected = kTRUE;
1262   if(evtCuts->IsEventSelectedRequired()) 
1263   { 
1264      // select events with at least 
1265      // one prompt track in acceptance
1266      // pT>0.5 GeV/c, |eta|<0.8 for the Cross Section studies
1267
1268      isEventSelected = AlidNdPtHelper::SelectEvent(esdEvent,esdTrackCuts);
1269      //printf("isEventSelected %d \n", isEventSelected);
1270   }
1271
1272   Bool_t isTrigAndVertex = isEventTriggered && isEventOK;
1273   Double_t vEventCount[3] = { static_cast<Double_t>(isEventTriggered), static_cast<Double_t>(isTrigAndVertex), static_cast<Double_t>(isEventSelected) };
1274   fEventCount->Fill(vEventCount);
1275
1276   //
1277   // correct event and track histograms
1278   //
1279   TObjArray *allChargedTracks=0;
1280   Int_t multRec=0, multRecTemp=0;
1281   Int_t *labelsRec=0;
1282   Bool_t isCosmic = kFALSE;
1283
1284
1285   if(isEventOK && isEventTriggered && isEventSelected)
1286   {
1287     // get all charged tracks
1288     allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());
1289     if(!allChargedTracks) return;
1290
1291     Int_t entries = allChargedTracks->GetEntries();
1292     labelsRec = new Int_t[entries];
1293
1294     // calculate mult of reconstructed tracks
1295
1296     for(Int_t i=0; i<entries;++i) 
1297     {
1298       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
1299       if(!track) continue;
1300       if(track->Charge()==0) continue;
1301
1302
1303       // only postive charged 
1304       if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0) 
1305         continue;
1306       
1307       // only negative charged 
1308       if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0) 
1309         continue;
1310
1311       // cosmics analysis
1312       isCosmic = kFALSE;
1313       if( GetParticleMode()==AlidNdPtHelper::kCosmic )
1314       {
1315           for(Int_t j=0; j<entries;++j) 
1316           {
1317             AliESDtrack *track1 = (AliESDtrack*)allChargedTracks->At(j);
1318             if(!track1) continue;
1319             if(track1->Charge()==0) continue;
1320
1321             if( esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track) && 
1322                 esdTrackCuts->AcceptTrack(track1) && accCuts->AcceptTrack(track1) ) 
1323             { 
1324               isCosmic = AlidNdPtHelper::IsCosmicTrack(track, track1);
1325             }
1326             if(isCosmic) 
1327             {
1328               Double_t vCosmicsHisto[3] = { track->Eta()+track1->Eta(), track->Phi()-track1->Phi(), track1->Pt() };
1329               fCosmicsHisto->Fill(vCosmicsHisto);
1330             }
1331           }
1332          
1333         if(!isCosmic) continue;
1334       }
1335
1336       if(esdTrackCuts->AcceptTrack(track)) 
1337       {
1338           if(accCuts->AcceptTrack(track)) multRecTemp++;
1339       }  
1340     }
1341
1342     //
1343     for(Int_t i=0; i<entries;++i) 
1344     {
1345       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);
1346       if(!track) continue;
1347       if(track->Charge()==0) continue;
1348
1349       // only postive charged 
1350       if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0) 
1351       continue;
1352       
1353       // only negative charged 
1354       if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0) 
1355       continue;
1356         
1357       // track-level corrections
1358       if(!esdTrackCuts->AcceptTrack(track))  continue;
1359       //if(GetAnalysisMode()==AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt && !AlidNdPtHelper::IsGoodImpPar(track)) continue;
1360
1361         // cosmics analysis
1362         isCosmic = kFALSE;
1363         if( GetParticleMode()==AlidNdPtHelper::kCosmic )
1364         {
1365           for(Int_t j=0; j<entries;++j) 
1366           {
1367             AliESDtrack *track1 = (AliESDtrack*)allChargedTracks->At(j);
1368             if(!track1) continue;
1369             if(track1->Charge()==0) continue;
1370
1371             if( esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track) && 
1372                 esdTrackCuts->AcceptTrack(track1) && accCuts->AcceptTrack(track1) ) 
1373             { 
1374               isCosmic = AlidNdPtHelper::IsCosmicTrack(track, track1);
1375             }
1376           }
1377           if(!isCosmic) continue;
1378         }
1379
1380         Bool_t isOK = kFALSE;
1381         Double_t x[3]; track->GetXYZ(x);
1382         Double_t b[3]; AliTracker::GetBxByBz(x,b);
1383
1384         //
1385         // if TPC-ITS hybrid tracking (kTPCITSHybrid)
1386         // replace track parameters with TPC-ony track parameters
1387         //
1388         if( GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybrid || GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtx || GetAnalysisMode() == AlidNdPtHelper::kTPCITSHybridTrackSPDvtxDCArPt) 
1389         {
1390           // Relate TPC-only tracks to Tracks or SPD vertex
1391           isOK = track->RelateToVertexTPCBxByBz(vtxESD, b, kVeryBig);
1392           if(!isOK) continue;
1393
1394           // replace esd track parameters with TPCinner
1395           AliExternalTrackParam  *tpcTrack  = new AliExternalTrackParam(*(track->GetTPCInnerParam()));
1396           if (!tpcTrack) return;
1397           track->Set(tpcTrack->GetX(),tpcTrack->GetAlpha(),tpcTrack->GetParameter(),tpcTrack->GetCovariance());
1398
1399           if(tpcTrack) delete tpcTrack; 
1400         } 
1401
1402         //
1403         if (GetAnalysisMode()==AlidNdPtHelper::kTPCSPDvtxUpdate || GetAnalysisMode()==AlidNdPtHelper::kTPCTrackSPDvtxUpdate) 
1404         {
1405            //
1406            // update track parameters
1407            //
1408            AliExternalTrackParam cParam;
1409            isOK = track->RelateToVertexTPC(vtxESD,esdEvent->GetMagneticField(),kVeryBig,&cParam);
1410            if(!isOK) continue;
1411            track->Set(cParam.GetX(),cParam.GetAlpha(),cParam.GetParameter(),cParam.GetCovariance());
1412
1413            if(accCuts->AcceptTrack(track)) { 
1414              FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxESD->GetZ(),multRecTemp,multMCTrueTracks); 
1415              labelsRec[multRec] = TMath::Abs(track->GetLabel());
1416              multRec++;
1417            }
1418          }
1419          else if (GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) 
1420          { 
1421            //
1422            // Replace rec with MC
1423            //
1424            if(accCuts->AcceptTrack(track)) { 
1425              FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxMC[2],multRecTemp, multMCTrueTracks); 
1426              labelsRec[multRec] = TMath::Abs(track->GetLabel());
1427              multRec++;
1428            }
1429          } 
1430          else  {
1431            //
1432            // all the rest tracking scenarios
1433            //
1434            if(accCuts->AcceptTrack(track)) { 
1435              FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxESD->GetZ(),multRecTemp, multMCTrueTracks); 
1436              labelsRec[multRec] = TMath::Abs(track->GetLabel());
1437              multRec++;
1438            }
1439          }
1440       }
1441
1442     // event-level corrections
1443     if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) 
1444     { 
1445       FillHistograms(AlidNdPtHelper::kRecEvents,vtxMC[2],multMBTracks);
1446     }
1447     else {
1448       FillHistograms(AlidNdPtHelper::kRecEvents,vtxESD->GetZ(),multMBTracks);
1449     }
1450
1451     // calculate meanPt from the event
1452     Double_t meanPtMult[8] = {0};  
1453     for (Int_t i = 0; i<8; i++) {
1454       if(!fCorrRecTrackMultHist1[i]) continue;
1455       TH1D *hp = (TH1D *)fCorrRecTrackPt1[i]->Projection(0);
1456       if(!hp) continue;
1457       meanPtMult[i] = hp->GetMean();    
1458       Double_t vCorrRecTrackMeanPtMultHist1[2] = {meanPtMult[i],static_cast<Double_t>(multRecTemp)};
1459       fCorrRecTrackMeanPtMultHist1[i]->Fill(vCorrRecTrackMeanPtMultHist1); 
1460       
1461       if( IsUseMCInfo() ) {
1462         Double_t vCorrRecTrackMeanPtTrueMultHist1[2] = {meanPtMult[i],static_cast<Double_t>(multMCTrueTracks)};
1463         fCorrRecTrackMeanPtTrueMultHist1[i]->Fill(vCorrRecTrackMeanPtTrueMultHist1); 
1464       }
1465
1466       // reset pt histo for the next event
1467       if(fCorrRecTrackPt1[i])  fCorrRecTrackPt1[i]->Reset();
1468       if(hp) delete hp;
1469     }
1470
1471     // control event histograms
1472     Double_t vRecEventHist1[3] = {vtxESD->GetX(),vtxESD->GetY(),vtxESD->GetZ()};
1473     fRecEventHist1->Fill(vRecEventHist1);
1474
1475     // correlation track multiplicity vs MB track multiplicity
1476     Double_t vRecEventMultHist1[3] = {static_cast<Double_t>(multRec), static_cast<Double_t>(multMBTracks)};
1477     fRecEventMultHist1->Fill(vRecEventMultHist1);
1478   }
1479
1480   // empty events corrections
1481   // no reconstructed zv
1482   if( isEventTriggered && multMBTracks==0 && isEventSelected ) 
1483   {
1484     if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) 
1485     {
1486       FillHistograms(AlidNdPtHelper::kTriggeredEvents,vtxMC[2],multMBTracks);
1487     }
1488     else {
1489       Double_t zv = fZvNorm->GetRandom();
1490       if(zv>evtCuts->GetMinZv() && zv<evtCuts->GetMaxZv())
1491         FillHistograms(AlidNdPtHelper::kTriggeredEvents,zv,multMBTracks);
1492     }
1493   }
1494
1495   if(IsUseMCInfo())  
1496   {
1497     if(!mcEvent) return; 
1498
1499     Bool_t isMCEventSelected = kTRUE;
1500     if(evtCuts->IsEventSelectedRequired()) 
1501     { 
1502       // select events with at least 
1503       // one MC primary track in acceptance
1504       // pT>0.5 GeV/c, |eta|<0.8 for the Cross Section studies
1505       isMCEventSelected = AlidNdPtHelper::SelectMCEvent(mcEvent);
1506       //printf("isMCEventSelected %d \n", isMCEventSelected);
1507     }
1508
1509     // select MC events 
1510     if(evtCuts->AcceptMCEvent(mcEvent) && isMCEventSelected)
1511     {
1512       //
1513       // event histograms
1514       //
1515       Double_t vMCEventMatrix[2] = {vtxMC[2],static_cast<Double_t>(multMBTracks)};
1516       fMCAllEventMultHist1->Fill(vMCEventMatrix);
1517
1518       if(evtType == AliPWG0Helper::kND) {
1519         fMCAllNDEventMultHist1->Fill(vMCEventMatrix);
1520       }
1521       if(evtType != AliPWG0Helper::kSD) {
1522         fMCAllNSDEventMultHist1->Fill(vMCEventMatrix);
1523       }
1524       if(isEventTriggered) {
1525         fMCTriggerMultHist1->Fill(vMCEventMatrix);
1526       }
1527       if(isEventTriggered && isEventOK) {
1528         fMCEventMultHist1->Fill(vMCEventMatrix);
1529       }
1530
1531       //
1532       // MC histograms for efficiency studies
1533       //
1534       Double_t sumPtMC = 0;
1535       Int_t nPart  = stack->GetNtrack();
1536       for (Int_t iMc = 0; iMc < nPart; ++iMc) 
1537       {
1538         // print MC stack info
1539         //AlidNdPtHelper::PrintMCInfo(stack,iMc);
1540
1541         TParticle* particle = stack->Particle(iMc);
1542         if (!particle)
1543         continue;
1544
1545         // only charged particles
1546           
1547         if(!particle->GetPDG()) continue;
1548         Double_t charge = particle->GetPDG()->Charge()/3.;
1549         if (TMath::Abs(charge) < 0.001)
1550           continue;
1551
1552         // only postive charged 
1553         if(GetParticleMode() == AlidNdPtHelper::kPlus && charge  < 0.) 
1554         continue;
1555       
1556         // only negative charged 
1557         if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.) 
1558         continue;
1559       
1560         // physical primary
1561         Bool_t prim = stack->IsPhysicalPrimary(iMc);
1562         if(!prim) continue;
1563
1564         // all primaries in acceptance
1565         if(!accCuts->AcceptTrack(particle)) continue;
1566
1567         Double_t gpt = particle->Pt();
1568         Double_t geta = particle->Eta();
1569
1570         // sum up pt in the event
1571         sumPtMC +=gpt; 
1572
1573         Double_t valMCAllTrackMultHist1[3] = {gpt,geta,static_cast<Double_t>(multRecTemp)};       
1574         Double_t valMCAllTrackTrueMultHist1[3] = {gpt,geta,static_cast<Double_t>(multMCTrueTracks)};      
1575         Double_t valMCAllTrackTrueMultHist2[3] = {gpt,static_cast<Double_t>(multRecTemp),static_cast<Double_t>(multMCTrueTracks)};        
1576
1577         fMCAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);
1578         fMCAllPrimTrackTrueMultHist1->Fill(valMCAllTrackTrueMultHist1);
1579         fMCAllPrimTrackTrueMultHist2->Fill(valMCAllTrackTrueMultHist2);
1580
1581         if(evtType == AliPWG0Helper::kND) {
1582           fMCNDEventAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);
1583           fMCNDEventAllPrimTrackTrueMultHist1->Fill(valMCAllTrackTrueMultHist1);
1584           fMCNDEventAllPrimTrackTrueMultHist2->Fill(valMCAllTrackTrueMultHist2);
1585         }
1586         if(evtType != AliPWG0Helper::kSD) {
1587           fMCNSDEventAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);
1588           fMCNSDEventAllPrimTrackTrueMultHist1->Fill(valMCAllTrackTrueMultHist1);
1589           fMCNSDEventAllPrimTrackTrueMultHist2->Fill(valMCAllTrackTrueMultHist2);
1590         }
1591         if(isEventTriggered) { 
1592           fMCTriggerPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);
1593           fMCTriggerPrimTrackTrueMultHist1->Fill(valMCAllTrackTrueMultHist1);
1594           fMCTriggerPrimTrackTrueMultHist2->Fill(valMCAllTrackTrueMultHist2);
1595         }
1596         if(isEventTriggered && isEventOK) { 
1597           fMCEventPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);
1598           fMCEventPrimTrackTrueMultHist1->Fill(valMCAllTrackTrueMultHist1);
1599           fMCEventPrimTrackTrueMultHist2->Fill(valMCAllTrackTrueMultHist2);
1600         }
1601       }
1602
1603       //
1604       // calculate <pt> in the event
1605       //
1606       Double_t meanPtMCMult = 0;
1607       Double_t meanPtMCTrueMult = 0;
1608       if(multRecTemp) { 
1609         meanPtMCMult = sumPtMC/multRecTemp; 
1610       }
1611       if(multMCTrueTracks) { 
1612         meanPtMCTrueMult = sumPtMC/multMCTrueTracks; 
1613       }
1614
1615       Double_t valMCMeanPtMult[2] = {static_cast<Double_t>(meanPtMCMult),static_cast<Double_t>(multRecTemp)};     
1616       Double_t valMCMeanPtTrueMult[2] = {static_cast<Double_t>(meanPtMCTrueMult),static_cast<Double_t>(multMCTrueTracks)};        
1617
1618       fMCAllPrimTrackMeanPtMult1->Fill(valMCMeanPtMult);
1619       fMCAllPrimTrackMeanPtTrueMult1->Fill(valMCMeanPtTrueMult);
1620
1621       if(evtType == AliPWG0Helper::kND) {
1622           fMCNDEventAllPrimTrackMeanPtMult1->Fill(valMCMeanPtMult);
1623           fMCNDEventAllPrimTrackMeanPtTrueMult1->Fill(valMCMeanPtTrueMult);
1624       }
1625       if(evtType != AliPWG0Helper::kSD) {
1626           fMCNSDEventAllPrimTrackMeanPtMult1->Fill(valMCMeanPtMult);
1627           fMCNSDEventAllPrimTrackMeanPtTrueMult1->Fill(valMCMeanPtTrueMult);
1628       }
1629       if(isEventTriggered) { 
1630           fMCTriggerPrimTrackMeanPtMult1->Fill(valMCMeanPtMult);
1631           fMCTriggerPrimTrackMeanPtTrueMult1->Fill(valMCMeanPtTrueMult);
1632       }
1633       if(isEventTriggered && isEventOK) { 
1634           fMCEventPrimTrackMeanPtMult1->Fill(valMCMeanPtMult);
1635           fMCEventPrimTrackMeanPtTrueMult1->Fill(valMCMeanPtTrueMult);
1636       }
1637     }
1638   } // end bUseMC
1639
1640   if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;
1641   if(labelsRec) delete [] labelsRec; labelsRec = 0;
1642
1643   if(!evtCuts->IsRecVertexRequired() && vtxESD != NULL) delete vtxESD;
1644 }
1645
1646 //_____________________________________________________________________________
1647 void AlidNdPtCorrection::FillHistograms(AlidNdPtHelper::EventObject eventObj, Double_t zv, Int_t multMBTracks) const
1648 {
1649   //
1650   // Fill corrected histograms
1651   //
1652
1653   Double_t vEventMatrix[2] = {zv,static_cast<Double_t>(multMBTracks)};
1654   //
1655   // Correct for efficiency 
1656   //
1657   if(eventObj == AlidNdPtHelper::kRecEvents && multMBTracks>0)  
1658   {
1659     Double_t corrToMBF = GetCorrFactZvMult(fCorrEventMatrix,zv,multMBTracks);
1660     Double_t corrToInelF = GetCorrFactZvMult(fCorrTriggerMBtoInelEventMatrix,zv,multMBTracks);
1661     Double_t corrToNDF = GetCorrFactZvMult(fCorrTriggerMBtoNDEventMatrix,zv,multMBTracks);
1662     Double_t corrToNSDF = GetCorrFactZvMult(fCorrTriggerMBtoNSDEventMatrix,zv,multMBTracks);
1663     //printf("corrToMBF %f, corrToInelF %f, corrToNDF %f corrToNSDF %f \n",corrToMBF,corrToInelF,corrToNDF,corrToNSDF);
1664
1665     fCorrRecEventHist1[0]->Fill(vEventMatrix);
1666     fCorrRecEventHist1[1]->Fill(vEventMatrix,corrToMBF);
1667     fCorrRecEventHist1[2]->Fill(vEventMatrix,corrToMBF*corrToInelF);
1668     fCorrRecEventHist1[3]->Fill(vEventMatrix,corrToMBF*corrToNDF);
1669     fCorrRecEventHist1[4]->Fill(vEventMatrix,corrToMBF*corrToNSDF);
1670   }
1671
1672   if(eventObj==AlidNdPtHelper::kTriggeredEvents && multMBTracks==0) // empty triggered events
1673   {
1674     Double_t factLHCBack = 1.;
1675     if(!IsUseMCInfo()) factLHCBack = fLHCBin0Background->GetBinContent(1); 
1676
1677
1678     Int_t bin = fZvEmptyEventsNorm->FindBin(zv);
1679     Double_t factZ = fZvEmptyEventsNorm->GetBinContent(bin);
1680
1681     Double_t corrToInelF0 = GetCorrFactZvMult(fCorrTriggerMBtoInelEventMatrix,zv,multMBTracks);
1682     Double_t corrToNDF0 = GetCorrFactZvMult(fCorrTriggerMBtoNDEventMatrix,zv,multMBTracks);
1683     Double_t corrToNSDF0 = GetCorrFactZvMult(fCorrTriggerMBtoNSDEventMatrix,zv,multMBTracks);
1684     //printf("factLHCBack %f, factZ %f, corrToInelF0 %f, corrToNDF0 %f, corrToNSDF0 %f \n",factLHCBack,factZ,corrToInelF0,corrToNDF0,corrToNSDF0);
1685
1686     fCorrRecEventHist2[0]->Fill(vEventMatrix);
1687     fCorrRecEventHist2[1]->Fill(vEventMatrix,factLHCBack*factZ);
1688     fCorrRecEventHist2[2]->Fill(vEventMatrix,factLHCBack*factZ*corrToInelF0);
1689     fCorrRecEventHist2[3]->Fill(vEventMatrix,factLHCBack*factZ*corrToNDF0);
1690     fCorrRecEventHist2[4]->Fill(vEventMatrix,factLHCBack*factZ*corrToNSDF0);
1691   }
1692 }
1693
1694 //_____________________________________________________________________________
1695 void AlidNdPtCorrection::FillHistograms(AliESDtrack * const esdTrack, AliStack * const stack, AlidNdPtHelper::TrackObject trackObj, Double_t zv, Int_t mult, Int_t trueMult) const
1696 {
1697   //
1698   // Fill ESD track and MC histograms 
1699   //
1700   if(!esdTrack) return;
1701
1702   //Float_t q = esdTrack->Charge();
1703   Float_t pt = esdTrack->Pt();
1704   Float_t eta = esdTrack->Eta();
1705   Float_t phi = esdTrack->Phi();
1706
1707   if(stack && GetAnalysisMode() == AlidNdPtHelper::kMCRec) 
1708   {
1709     Int_t label = TMath::Abs(esdTrack->GetLabel());
1710    
1711     TParticle* particle = stack->Particle(label);
1712     if(!particle) return;
1713    
1714     if(!particle->GetPDG()) return;
1715     Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3
1716     if(TMath::Abs(gq)<0.001) return;
1717     Float_t gpt = particle->Pt();
1718     Float_t geta = particle->Eta();
1719     Float_t gphi = particle->Phi();
1720
1721     // replace reconstructed values with MC
1722     pt = gpt;
1723     eta = geta;
1724     phi = gphi;
1725   }
1726
1727   //
1728   // Fill histograms
1729   //
1730   Double_t values[3] = {pt,eta,phi};      
1731   fRecTrackHist1[trackObj]->Fill(values);
1732
1733   //
1734   // Correct for contamination and efficiency 
1735   //
1736   if(trackObj == AlidNdPtHelper::kRecTracks || GetAnalysisMode() == AlidNdPtHelper::kMCRec)  
1737   {
1738     // track level corrections
1739     Double_t trackEffF = 1.0;  
1740     if(pt < 2.6) trackEffF = GetCorrFactZvPtEta(fCorrTrackMatrix,zv,pt,eta);
1741     else trackEffF = GetCorrFactZvPtEta(fCorrHighPtTrackMatrix,zv,pt,eta);
1742
1743     Double_t trackContF = GetContFactZvPtEta(fContTrackMatrix,zv,pt,eta);
1744     Double_t multTrackContF = GetContFactZvPtEta(fContMultTrackMatrix,zv,pt,eta);
1745     //printf("zv %f, pt %f, eta %f \n",zv,pt,eta);
1746     //printf("trackEffF %f, trackContF %f, multTrackContF %f \n", trackEffF, trackContF, multTrackContF);
1747    
1748     // track-event level corrections
1749     Double_t vertexEffF = GetCorrFactZvPtEta(fCorrTrackEventMatrix,zv,pt,eta);
1750     Double_t trigMBToInel = GetCorrFactZvPtEta(fCorrTriggerMBtoInelTrackEventMatrix,zv,pt,eta);  
1751     Double_t trigMBToND = GetCorrFactZvPtEta(fCorrTriggerMBtoNDTrackEventMatrix,zv,pt,eta);
1752     Double_t trigMBToNSD = GetCorrFactZvPtEta(fCorrTriggerMBtoNSDTrackEventMatrix,zv,pt,eta);
1753     //printf("vertexEffF %f, trigMBToInel %f, trigMBToNSD %f \n", vertexEffF, trigMBToInel, trigMBToNSD);
1754     
1755     Double_t corrF[8] = { 1.0, 
1756                           trackContF,
1757                           trackContF*trackEffF,
1758                           trackContF*trackEffF*multTrackContF,
1759                           trackContF*trackEffF*multTrackContF*vertexEffF,
1760                           trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToInel,
1761                           trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToND,
1762                           trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToNSD
1763                          }; 
1764  
1765     // Fill histograms
1766     Double_t valCorrRecTrackMultHist1[3] = {pt,eta,static_cast<Double_t>(mult)};          
1767     Double_t valCorrRecTrackPt1[1] = {pt};        
1768     for(Int_t i=0; i<8; i++) {
1769       fCorrRecTrackMultHist1[i]->Fill(valCorrRecTrackMultHist1,corrF[i]);
1770       fCorrRecTrackPt1[i]->Fill(valCorrRecTrackPt1,corrF[i]);
1771
1772       if( IsUseMCInfo() ) {
1773         Double_t valCorrRecTrackTrueMultHist1[3] = {pt,eta,static_cast<Double_t>(trueMult)};      
1774         Double_t valCorrRecTrackTrueMultHist2[3] = {pt,static_cast<Double_t>(mult),static_cast<Double_t>(trueMult)};      
1775
1776         fCorrRecTrackTrueMultHist1[i]->Fill(valCorrRecTrackTrueMultHist1,corrF[i]);
1777         fCorrRecTrackTrueMultHist2[i]->Fill(valCorrRecTrackTrueMultHist2,corrF[i]);
1778       } 
1779     }
1780   }
1781 }
1782
1783 void AlidNdPtCorrection::FillHistograms(AliStack * const stack, Int_t /*label*/, AlidNdPtHelper::TrackObject /*trackObj*/, Int_t /*mult*/) const
1784 {
1785   // Fill MC histograms
1786   if(!stack) return;
1787
1788   /*
1789   TParticle* particle = stack->Particle(label);
1790   if(!particle) return;
1791
1792   Int_t mother_pdg = -1;
1793   TParticle* mother = 0;
1794
1795   //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);
1796   Int_t motherLabel = particle->GetMother(0); 
1797   if(motherLabel>0) mother = stack->Particle(motherLabel);
1798   if(mother) mother_pdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only
1799   Int_t mech = particle->GetUniqueID(); // production mechanism
1800
1801   if(!particle->GetPDG()) return;
1802   Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 
1803   Float_t gpt = particle->Pt();
1804   Float_t qgpt = particle->Pt() * gq;
1805   Float_t geta = particle->Eta();
1806   Float_t gphi = particle->Phi();
1807   Float_t gpz = particle->Pz();
1808
1809   Bool_t prim = stack->IsPhysicalPrimary(label);
1810   Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();
1811
1812   Int_t pid=-1;
1813   if (TMath::Abs(particle->GetPdgCode()) == kElectron)         { pid = 0; }
1814     else if (TMath::Abs(particle->GetPdgCode()) == kMuonMinus) { pid = 1; }
1815     else if (TMath::Abs(particle->GetPdgCode()) == kPiPlus)    { pid = 2; }
1816     else if (TMath::Abs(particle->GetPdgCode()) == kKPlus)     { pid = 3; }
1817     else if (TMath::Abs(particle->GetPdgCode()) == kProton)    { pid = 4; }
1818     else                                                       { pid = 5; }
1819     */
1820
1821   //if(!prim) printf("prim_mother %d, mother %d, particle %d, production mech %d\n",prim_mother->GetPdgCode(),mother->GetPdgCode(), particle->GetPdgCode(),mech);
1822   
1823 }
1824
1825 //_____________________________________________________________________________
1826 Double_t AlidNdPtCorrection::GetCorrFactZvPtEta(THnSparse * const hist, Double_t zv, Double_t pt, Double_t eta) const {
1827 // return correction factor F(zv,pt,eta)
1828
1829  if(!hist) return 1.;
1830
1831  //
1832  TAxis *ax = hist->GetAxis(0);
1833  TAxis *ay = hist->GetAxis(1);
1834  TAxis *az = hist->GetAxis(2);
1835
1836  Int_t binx = ax->FindBin(zv);
1837  Int_t biny = ay->FindBin(pt);
1838  Int_t binz = az->FindBin(eta);
1839  Int_t dim[3] = {binx,biny,binz};
1840
1841  Double_t fact  = hist->GetBinContent(dim);  
1842
1843 return fact;
1844 }
1845
1846 //_____________________________________________________________________________
1847 Double_t AlidNdPtCorrection::GetContFactZvPtEta(THnSparse * const hist, Double_t zv, Double_t pt, Double_t eta) const {
1848 // return contamination correction factor F(zv,pt,eta)
1849
1850  if(!hist) return 1.0;
1851
1852  //
1853  TAxis *ax = hist->GetAxis(0);
1854  TAxis *ay = hist->GetAxis(1);
1855  TAxis *az = hist->GetAxis(2);
1856
1857  Int_t binx = ax->FindBin(zv);
1858  Int_t biny = ay->FindBin(pt);
1859  Int_t binz = az->FindBin(eta);
1860  Int_t dim[3] = {binx,biny,binz};
1861
1862  //
1863  //  additional correction for secondary 
1864  //  particles with strangeness (data driven)
1865  // 
1866  Double_t corrFact = 1.;
1867  if(!IsUseMCInfo()) corrFact = AlidNdPtHelper::GetStrangenessCorrFactor(pt);
1868  //printf("pt %f, corrFact %f \n", pt, corrFact);
1869
1870  Double_t fact  = 1.0 - corrFact*hist->GetBinContent(dim);  
1871  //Double_t fact  = hist->GetBinContent(dim);  
1872
1873 return fact;
1874 }
1875
1876 //_____________________________________________________________________________
1877 Double_t AlidNdPtCorrection::GetCorrFactZvMult(THnSparse * const hist, Double_t zv, Int_t mult) const {
1878 // return correction factor F(zv,mult)
1879
1880  if(!hist) return 1.;
1881
1882  TAxis *ax = hist->GetAxis(0);
1883  TAxis *ay = hist->GetAxis(1);
1884  Int_t binx = ax->FindBin(zv);
1885  Int_t biny = ay->FindBin(mult);
1886  Int_t dim[2] = {binx,biny};
1887
1888  Double_t fact  = hist->GetBinContent(dim);  
1889
1890
1891 return fact;
1892 }
1893
1894 //_____________________________________________________________________________
1895 Double_t AlidNdPtCorrection::GetContFactZvMult(THnSparse * const hist, Double_t zv, Int_t mult) const {
1896 // return contamination correction factor F(zv,mult)
1897
1898  if(!hist) return 1.;
1899
1900  TAxis *ax = hist->GetAxis(0);
1901  TAxis *ay = hist->GetAxis(1);
1902  Int_t binx = ax->FindBin(zv);
1903  Int_t biny = ay->FindBin(mult);
1904  Int_t dim[2] = {binx,biny};
1905  Double_t fact  = 1.0-hist->GetBinContent(dim);  
1906
1907 return fact;
1908 }
1909
1910 //_____________________________________________________________________________
1911 Long64_t AlidNdPtCorrection::Merge(TCollection* const list) 
1912 {
1913   // Merge list of objects (needed by PROOF)
1914
1915   if (!list)
1916   return 0;
1917
1918   if (list->IsEmpty())
1919   return 1;
1920
1921   TIterator* iter = list->MakeIterator();
1922   TObject* obj = 0;
1923
1924   // collection of generated histograms
1925
1926   // physics selection
1927   TList *collPhysSelection = new TList;
1928
1929   Int_t count=0;
1930   while((obj = iter->Next()) != 0) {
1931     AlidNdPtCorrection* entry = dynamic_cast<AlidNdPtCorrection*>(obj);
1932     if (entry == 0) continue; 
1933   
1934     collPhysSelection->Add(entry->GetPhysicsTriggerSelection());
1935
1936     fEventCount->Add(entry->fEventCount);
1937
1938     fMCEventHist1->Add(entry->fMCEventHist1);
1939     fRecEventHist1->Add(entry->fRecEventHist1);
1940     fRecEventMultHist1->Add(entry->fRecEventMultHist1);
1941
1942     fMCAllEventMultHist1->Add(entry->fMCAllEventMultHist1);
1943     fMCAllNDEventMultHist1->Add(entry->fMCAllNDEventMultHist1);
1944     fMCAllNSDEventMultHist1->Add(entry->fMCAllNSDEventMultHist1);
1945     fMCTriggerMultHist1->Add(entry->fMCTriggerMultHist1);
1946     fMCEventMultHist1->Add(entry->fMCEventMultHist1);
1947
1948     fMCAllPrimTrackMultHist1->Add(entry->fMCAllPrimTrackMultHist1);
1949     fMCNDEventAllPrimTrackMultHist1->Add(entry->fMCNDEventAllPrimTrackMultHist1);
1950     fMCNSDEventAllPrimTrackMultHist1->Add(entry->fMCNSDEventAllPrimTrackMultHist1);
1951     fMCTriggerPrimTrackMultHist1->Add(entry->fMCTriggerPrimTrackMultHist1);
1952     fMCEventPrimTrackMultHist1->Add(entry->fMCEventPrimTrackMultHist1);
1953
1954     fMCAllPrimTrackTrueMultHist1->Add(entry->fMCAllPrimTrackTrueMultHist1);
1955     fMCNDEventAllPrimTrackTrueMultHist1->Add(entry->fMCNDEventAllPrimTrackTrueMultHist1);
1956     fMCNSDEventAllPrimTrackTrueMultHist1->Add(entry->fMCNSDEventAllPrimTrackTrueMultHist1);
1957     fMCTriggerPrimTrackTrueMultHist1->Add(entry->fMCTriggerPrimTrackTrueMultHist1);
1958     fMCEventPrimTrackTrueMultHist1->Add(entry->fMCEventPrimTrackTrueMultHist1);
1959
1960     fMCAllPrimTrackTrueMultHist2->Add(entry->fMCAllPrimTrackTrueMultHist2);
1961     fMCNDEventAllPrimTrackTrueMultHist2->Add(entry->fMCNDEventAllPrimTrackTrueMultHist2);
1962     fMCNSDEventAllPrimTrackTrueMultHist2->Add(entry->fMCNSDEventAllPrimTrackTrueMultHist2);
1963     fMCTriggerPrimTrackTrueMultHist2->Add(entry->fMCTriggerPrimTrackTrueMultHist2);
1964     fMCEventPrimTrackTrueMultHist2->Add(entry->fMCEventPrimTrackTrueMultHist2);
1965
1966     fMCAllPrimTrackMeanPtMult1->Add(entry->fMCAllPrimTrackMeanPtMult1);
1967     fMCNDEventAllPrimTrackMeanPtMult1->Add(entry->fMCNDEventAllPrimTrackMeanPtMult1);
1968     fMCNSDEventAllPrimTrackMeanPtMult1->Add(entry->fMCNSDEventAllPrimTrackMeanPtMult1);
1969     fMCTriggerPrimTrackMeanPtMult1->Add(entry->fMCTriggerPrimTrackMeanPtMult1);
1970     fMCEventPrimTrackMeanPtMult1->Add(entry->fMCEventPrimTrackMeanPtMult1);
1971
1972     fMCAllPrimTrackMeanPtTrueMult1->Add(entry->fMCAllPrimTrackMeanPtTrueMult1);
1973     fMCNDEventAllPrimTrackMeanPtTrueMult1->Add(entry->fMCNDEventAllPrimTrackMeanPtTrueMult1);
1974     fMCNSDEventAllPrimTrackMeanPtTrueMult1->Add(entry->fMCNSDEventAllPrimTrackMeanPtTrueMult1);
1975     fMCTriggerPrimTrackMeanPtTrueMult1->Add(entry->fMCTriggerPrimTrackMeanPtTrueMult1);
1976     fMCEventPrimTrackMeanPtTrueMult1->Add(entry->fMCEventPrimTrackMeanPtTrueMult1);
1977
1978     fCosmicsHisto->Add(entry->fCosmicsHisto);
1979
1980     for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {
1981       fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);
1982     }
1983
1984     for(Int_t i=0; i<8; i++) {
1985       fCorrRecTrackMultHist1[i]->Add(entry->fCorrRecTrackMultHist1[i]);
1986       fCorrRecTrackTrueMultHist1[i]->Add(entry->fCorrRecTrackTrueMultHist1[i]);
1987       fCorrRecTrackTrueMultHist2[i]->Add(entry->fCorrRecTrackTrueMultHist2[i]);
1988
1989       fCorrRecTrackMeanPtMultHist1[i]->Add(entry->fCorrRecTrackMeanPtMultHist1[i]);
1990       fCorrRecTrackMeanPtTrueMultHist1[i]->Add(entry->fCorrRecTrackMeanPtTrueMultHist1[i]);
1991
1992       fCorrRecTrackPt1[i]->Add(entry->fCorrRecTrackPt1[i]);
1993     }
1994
1995     for(Int_t i=0; i<5; i++) {
1996       fCorrRecEventHist1[i]->Add(entry->fCorrRecEventHist1[i]);
1997       fCorrRecEventHist2[i]->Add(entry->fCorrRecEventHist2[i]);
1998     }
1999
2000   count++;
2001   }
2002
2003   //
2004   AliPhysicsSelection *trigSelection = GetPhysicsTriggerSelection();
2005   trigSelection->Merge(collPhysSelection);
2006   if(collPhysSelection) delete collPhysSelection;
2007
2008 return count;
2009 }
2010  
2011 //____________________________________________________________________________
2012 Int_t AlidNdPtCorrection::GetTrueMult(THnSparse * const hist, Int_t mult) const
2013 {
2014 //
2015 // get multiplicity of primary particles
2016 //
2017  if(!hist) return 0;
2018  Int_t trueMult = 0;
2019
2020  // 0 bins exluded
2021  TAxis *ax = hist->GetAxis(0);
2022  TAxis *ay = hist->GetAxis(1);
2023  ax->SetRange(1,ax->GetNbins());
2024  ay->SetRange(1,ay->GetNbins());
2025
2026  // measured mult
2027  ax->SetRangeUser((Float_t)mult,(Float_t)mult); 
2028
2029  // get true multiplicity
2030  TH1D *h1 = (TH1D *)hist->Projection(1);
2031  trueMult = (Int_t)h1->GetMean();
2032
2033  return trueMult;
2034 }
2035
2036 //_____________________________________________________________________________
2037 void AlidNdPtCorrection::Analyse() 
2038 {
2039   // Analyse histograms
2040   //
2041   TH1::AddDirectory(kFALSE);
2042   TH1 *h = 0, *hs=0, *hsc=0; 
2043   TH2 *h2D = 0; 
2044
2045   TObjArray *aFolderObj = new TObjArray;
2046   if(!aFolderObj) return;
2047
2048   //
2049   // get cuts
2050   //
2051   AlidNdPtEventCuts *evtCuts = GetEventCuts(); 
2052   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); 
2053   AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); 
2054
2055   if(!evtCuts || !accCuts || !esdTrackCuts) {
2056     Error("AlidNdPtCutAnalysis::Analyse()", "cuts not available");
2057     return;
2058   }
2059
2060   //
2061   // set min and max values
2062   //
2063   //Double_t minPt = accCuts->GetMinPt();
2064   //Double_t maxPt = accCuts->GetMaxPt();
2065   Double_t minEta = accCuts->GetMinEta();
2066   Double_t maxEta = accCuts->GetMaxEta()-0.00001;
2067  
2068   printf("minEta %f, maxEta %f \n",minEta, maxEta);
2069
2070   //
2071   // LHC backgraund in all and 0-bins
2072   //
2073   AliPhysicsSelection *trigSel = GetPhysicsTriggerSelection();
2074   trigSel->SaveHistograms("physics_selection");
2075
2076   //
2077   // cosmics background histo
2078   //
2079   h2D = fCosmicsHisto->Projection(0,1);
2080   h2D->SetName("deta_vs_dphi_cosmics");
2081   aFolderObj->Add(h2D);
2082
2083   //
2084   // event level 
2085   //
2086   h = fCorrRecEventHist1[0]->Projection(1);
2087   h->SetName("mult_event_not_corrected");
2088   aFolderObj->Add(h);
2089
2090   h = fCorrRecEventHist1[1]->Projection(1);
2091   h->SetName("mult_event_vertex_corrected");
2092   aFolderObj->Add(h);
2093
2094   h = fCorrRecEventHist1[2]->Projection(1);
2095   h->SetName("mult_trigger_vertex_corrected");
2096   aFolderObj->Add(h);
2097
2098   h = fCorrRecEventHist1[3]->Projection(1);
2099   h->SetName("mult_ND_trigger_vertex_corrected");
2100   aFolderObj->Add(h);
2101
2102   h = fCorrRecEventHist1[4]->Projection(1);
2103   h->SetName("mult_NSD_trigger_vertex_corrected");
2104   aFolderObj->Add(h);
2105
2106   // empty events
2107   h = fCorrRecEventHist2[0]->Projection(1);
2108   h->SetName("mult_empty_event_not_corrected");
2109   aFolderObj->Add(h);
2110
2111   h = fCorrRecEventHist2[1]->Projection(1);
2112   h->SetName("mult_empty_event_vertex_corrected");
2113   aFolderObj->Add(h);
2114
2115   h = fCorrRecEventHist2[2]->Projection(1);
2116   h->SetName("mult_empty_trigger_vertex_corrected");
2117   aFolderObj->Add(h);
2118
2119   h = fCorrRecEventHist2[3]->Projection(1);
2120   h->SetName("mult_empty_ND_trigger_vertex_corrected");
2121   aFolderObj->Add(h);
2122
2123   h = fCorrRecEventHist2[4]->Projection(1);
2124   h->SetName("mult_empty_NSD_trigger_vertex_corrected");
2125   aFolderObj->Add(h);
2126  
2127   //
2128   // MC available
2129   //
2130   if(IsUseMCInfo()) {
2131
2132   // mc 
2133   h = fMCAllEventMultHist1->Projection(1);
2134   h->SetName("mc_mult_event_acc_prim");
2135   aFolderObj->Add(h);
2136
2137   h = fMCAllNDEventMultHist1->Projection(1);
2138   h->SetName("mc_mult_ND_event_acc_prim");
2139   aFolderObj->Add(h);
2140
2141   h = fMCAllNSDEventMultHist1->Projection(1);
2142   h->SetName("mc_mult_NSD_event_acc_prim");
2143   aFolderObj->Add(h);
2144
2145   h = fMCTriggerMultHist1->Projection(1);
2146   h->SetName("mc_mult_trigger_acc_prim");
2147   aFolderObj->Add(h);
2148
2149   h = fMCEventMultHist1->Projection(1);
2150   h->SetName("mc_mult_trigger_event_acc_prim");
2151   aFolderObj->Add(h);
2152
2153
2154   //
2155   // track level
2156   //
2157   
2158   // limit eta range
2159   for(Int_t i=0;i<8;i++) { 
2160       //fCorrRecTrackMultHist1[i]->GetAxis(0)->SetRangeUser(minPt,maxPt);
2161       //fCorrRecTrackMultHist1[i]->GetAxis(1)->SetRangeUser(minEta,maxEta);
2162   }
2163   //fMCAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);
2164   //fMCAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);
2165
2166   //fMCNDEventAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);
2167   //fMCNDEventAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);
2168
2169   //fMCNSDEventAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);
2170   //fMCNSDEventAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);
2171
2172   //fMCTriggerPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);
2173   //fMCTriggerPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);
2174
2175   //fMCEventPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);
2176   //fMCEventPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);
2177
2178   } // end use MC info 
2179   
2180   //
2181   h2D = fCorrRecTrackMultHist1[3]->Projection(1,0);
2182   h2D->SetName("pt_eta_rec_track_mult_eff_cont_corrected");
2183   aFolderObj->Add(h2D);
2184
2185   h2D = fCorrRecTrackMultHist1[4]->Projection(1,0);
2186   h2D->SetName("pt_eta_rec_event_track_mult_eff_cont_corrected");
2187   aFolderObj->Add(h2D);
2188
2189   h2D = fCorrRecTrackMultHist1[5]->Projection(1,0);
2190   h2D->SetName("pt_eta_rec_trig_event_track_mult_eff_cont_corrected");
2191   aFolderObj->Add(h2D);
2192
2193   h2D = fCorrRecTrackMultHist1[6]->Projection(1,0);
2194   h2D->SetName("pt_eta_rec_ND_trig_event_track_mult_eff_cont_corrected");
2195   aFolderObj->Add(h2D);
2196
2197   h2D = fCorrRecTrackMultHist1[7]->Projection(1,0);
2198   h2D->SetName("pt_eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");
2199   aFolderObj->Add(h2D);
2200
2201
2202   //
2203   h2D = fCorrRecTrackMultHist1[3]->Projection(2,0);
2204   h2D->SetName("pt_mult_rec_track_mult_eff_cont_corrected");
2205   aFolderObj->Add(h2D);
2206
2207   h2D = fCorrRecTrackMultHist1[4]->Projection(2,0);
2208   h2D->SetName("pt_mult_rec_event_track_mult_eff_cont_corrected");
2209   aFolderObj->Add(h2D);
2210
2211   h2D = fCorrRecTrackMultHist1[5]->Projection(2,0);
2212   h2D->SetName("pt_mult_rec_trig_event_track_mult_eff_cont_corrected");
2213   aFolderObj->Add(h2D);
2214
2215   h2D = fCorrRecTrackMultHist1[6]->Projection(2,0);
2216   h2D->SetName("pt_mult_rec_ND_trig_event_track_mult_eff_cont_corrected");
2217   aFolderObj->Add(h2D);
2218
2219   h2D = fCorrRecTrackMultHist1[7]->Projection(2,0);
2220   h2D->SetName("pt_mult_rec_NSD_trig_event_track_mult_eff_cont_corrected");
2221   aFolderObj->Add(h2D);
2222
2223   // pt axis
2224
2225   h = fCorrRecTrackMultHist1[0]->Projection(0);
2226   h->SetName("pt_rec_track_not_corrected");
2227   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2228   hs->SetName("pt_rec_track_not_corrected_s");
2229   aFolderObj->Add(hs);
2230
2231   //
2232   h = fCorrRecTrackMultHist1[1]->Projection(0);
2233   h->SetName("pt_rec_track_cont_corrected");
2234   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2235   hs->SetName("pt_rec_track_cont_corrected_s");
2236   aFolderObj->Add(hs);
2237
2238   hsc = (TH1D*)hs->Clone("pt_rec_track_cont_corr_fact");
2239   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_not_corrected_s"));
2240   aFolderObj->Add(hsc);
2241
2242   //
2243   h = fCorrRecTrackMultHist1[2]->Projection(0);
2244   h->SetName("pt_rec_track_eff_cont_corrected");
2245   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2246   hs->SetName("pt_rec_track_eff_cont_corrected_s");
2247   aFolderObj->Add(hs);
2248
2249   hsc = (TH1D*)hs->Clone("pt_rec_track_eff_corr_fact");
2250   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_cont_corrected_s"));
2251   aFolderObj->Add(hsc);
2252
2253   //
2254   h = fCorrRecTrackMultHist1[3]->Projection(0);
2255   h->SetName("pt_rec_track_mult_eff_cont_corrected");
2256   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2257   hs->SetName("pt_rec_track_mult_eff_cont_corrected_s");
2258   aFolderObj->Add(hs);
2259
2260   hsc = (TH1D*)hs->Clone("pt_rec_track_mult_corr_fact");
2261   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_eff_cont_corrected_s"));
2262   aFolderObj->Add(hsc);
2263
2264   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2265   hsc = (TH1D*)hs->Clone("pt_rec_track_mult_eff_cont_corr_fact");
2266   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_not_corrected_s"));
2267   aFolderObj->Add(hsc);
2268
2269   hsc = (TH1D*)hs->Clone();
2270   hsc->SetName("pt_rec_track_mult_eff_cont_corrected_s_norm");
2271   hsc->Scale(1./(fCorrRecEventHist1[0]->Projection(1)->Integral()));
2272   aFolderObj->Add(hsc);
2273
2274   //
2275   h = fCorrRecTrackMultHist1[4]->Projection(0);
2276   h->SetName("pt_rec_event_track_mult_eff_cont_corrected");
2277   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2278   hs->SetName("pt_rec_event_track_mult_eff_cont_corrected_s");
2279   aFolderObj->Add(hs);
2280
2281   hsc = (TH1D*)hs->Clone();
2282   hsc->SetName("pt_rec_event_track_mult_eff_cont_corrected_s_norm");
2283   hsc->Scale(1./(fCorrRecEventHist1[1]->Projection(1)->Integral()+fCorrRecEventHist2[1]->Projection(1)->Integral()));
2284   aFolderObj->Add(hsc);
2285
2286   //
2287   h = fCorrRecTrackMultHist1[5]->Projection(0);
2288   h->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected");
2289   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2290   hs->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s");
2291   aFolderObj->Add(hs);
2292
2293   hsc = (TH1D*)hs->Clone();
2294   hsc->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm");
2295   hsc->Scale(1./(fCorrRecEventHist1[2]->Projection(1)->Integral() + fCorrRecEventHist2[2]->Projection(1)->Integral()));
2296   aFolderObj->Add(hsc);
2297
2298    // positive eta
2299   fCorrRecTrackMultHist1[5]->GetAxis(1)->SetRangeUser(0., maxEta);
2300
2301   h = fCorrRecTrackMultHist1[5]->Projection(0);
2302   h->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_posEta");
2303   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2304   hs->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_posEta");
2305   aFolderObj->Add(hs);
2306
2307   hsc = (TH1D*)hs->Clone();
2308   hsc->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm_posEta");
2309   hsc->Scale(1./(fCorrRecEventHist1[2]->Projection(1)->Integral()+fCorrRecEventHist2[2]->Projection(1)->Integral()));
2310   aFolderObj->Add(hsc);
2311
2312   // negative eta
2313   fCorrRecTrackMultHist1[5]->GetAxis(1)->SetRangeUser(minEta, -0.00001);
2314
2315   h = fCorrRecTrackMultHist1[5]->Projection(0);
2316   h->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_negEta");
2317   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2318   hs->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_negEta");
2319   aFolderObj->Add(hs);
2320
2321   hsc = (TH1D*)hs->Clone();
2322   hsc->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm_negEta");
2323   hsc->Scale(1./(fCorrRecEventHist1[2]->Projection(1)->Integral()+fCorrRecEventHist2[2]->Projection(1)->Integral()));
2324   aFolderObj->Add(hsc);
2325
2326   fCorrRecTrackMultHist1[5]->GetAxis(1)->SetRange(1, fCorrRecTrackMultHist1[5]->GetAxis(1)->GetNbins());
2327
2328   //
2329   h = fCorrRecTrackMultHist1[6]->Projection(0);
2330   h->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected");
2331   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2332   hs->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s");
2333   aFolderObj->Add(hs);
2334
2335   hsc = (TH1D*)hs->Clone();
2336   hsc->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s_norm");
2337   hsc->Scale(1./(fCorrRecEventHist1[3]->Projection(1)->Integral()+fCorrRecEventHist2[3]->Projection(1)->Integral()));
2338   aFolderObj->Add(hsc);
2339
2340   //
2341   h = fCorrRecTrackMultHist1[7]->Projection(0);
2342   h->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected");
2343   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2344   hs->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s");
2345   aFolderObj->Add(hs);
2346
2347   hsc = (TH1D*)hs->Clone();
2348   hsc->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm");
2349   hsc->Scale(1./(fCorrRecEventHist1[4]->Projection(1)->Integral() + fCorrRecEventHist2[4]->Projection(1)->Integral()));
2350   aFolderObj->Add(hsc);
2351
2352   //
2353   // positive eta
2354   //
2355   fCorrRecTrackMultHist1[7]->GetAxis(1)->SetRangeUser(0., maxEta);
2356
2357   h = fCorrRecTrackMultHist1[7]->Projection(0);
2358   h->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_posEta");
2359   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2360   hs->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_posEta");
2361   aFolderObj->Add(hs);
2362
2363   hsc = (TH1D*)hs->Clone();
2364   hsc->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm_posEta");
2365   hsc->Scale(1./(fCorrRecEventHist1[4]->Projection(1)->Integral()+fCorrRecEventHist2[4]->Projection(1)->Integral()));
2366   aFolderObj->Add(hsc);
2367
2368   //
2369   // negative eta
2370   //
2371   fCorrRecTrackMultHist1[7]->GetAxis(1)->SetRangeUser(minEta, -0.00001);
2372
2373   h = fCorrRecTrackMultHist1[7]->Projection(0);
2374   h->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_negEta");
2375   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2376   hs->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_negEta");
2377   aFolderObj->Add(hs);
2378
2379   hsc = (TH1D*)hs->Clone();
2380   hsc->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm_negEta");
2381   hsc->Scale(1./(fCorrRecEventHist1[4]->Projection(1)->Integral()+fCorrRecEventHist2[4]->Projection(1)->Integral()));
2382   aFolderObj->Add(hsc);
2383
2384   fCorrRecTrackMultHist1[7]->GetAxis(1)->SetRange(1, fCorrRecTrackMultHist1[7]->GetAxis(1)->GetNbins());
2385
2386   // eta axis
2387   h = fCorrRecTrackMultHist1[0]->Projection(1);
2388   h->SetName("eta_rec_track_not_corrected");
2389   aFolderObj->Add(h);
2390   
2391   h = fCorrRecTrackMultHist1[1]->Projection(1);
2392   h->SetName("eta_rec_track_cont_corrected");
2393   aFolderObj->Add(h);
2394
2395   h = fCorrRecTrackMultHist1[2]->Projection(1);
2396   h->SetName("eta_rec_track_eff_cont_corrected");
2397   aFolderObj->Add(h);
2398
2399   h = fCorrRecTrackMultHist1[3]->Projection(1);
2400   h->SetName("eta_rec_track_mult_eff_cont_corrected");
2401   aFolderObj->Add(h);
2402
2403   h = fCorrRecTrackMultHist1[4]->Projection(1);
2404   h->SetName("eta_rec_event_track_mult_eff_cont_corrected");
2405   aFolderObj->Add(h);
2406
2407   h = fCorrRecTrackMultHist1[5]->Projection(1);
2408   h->SetName("eta_rec_trig_event_track_mult_eff_cont_corrected");
2409   aFolderObj->Add(h);
2410
2411   h = fCorrRecTrackMultHist1[6]->Projection(1);
2412   h->SetName("eta_rec_ND_trig_event_track_mult_eff_cont_corrected");
2413   aFolderObj->Add(h);
2414
2415   h = fCorrRecTrackMultHist1[7]->Projection(1);
2416   h->SetName("eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");
2417   aFolderObj->Add(h);
2418
2419
2420   //
2421   // MC available
2422   //
2423   if(IsUseMCInfo()) {
2424
2425   //
2426   h2D = fMCAllPrimTrackMultHist1->Projection(2,0);
2427   h2D->SetName("mc_all_pt_mult_acc_prim");
2428   aFolderObj->Add(h2D);
2429
2430   h2D = fMCAllPrimTrackMultHist1->Projection(1,0);
2431   h2D->SetName("mc_all_eta_pt_acc_prim");
2432   aFolderObj->Add(h2D);
2433
2434   h2D = fMCNDEventAllPrimTrackMultHist1->Projection(2,0);
2435   h2D->SetName("mc_ND_all_pt_mult_acc_prim");
2436   aFolderObj->Add(h2D);
2437
2438   h2D = fMCNDEventAllPrimTrackMultHist1->Projection(1,0);
2439   h2D->SetName("mc_ND_all_eta_pt_acc_prim");
2440   aFolderObj->Add(h2D);
2441
2442   h2D = fMCNSDEventAllPrimTrackMultHist1->Projection(2,0);
2443   h2D->SetName("mc_NSD_all_pt_mult_acc_prim");
2444   aFolderObj->Add(h2D);
2445
2446   h2D = fMCNSDEventAllPrimTrackMultHist1->Projection(1,0);
2447   h2D->SetName("mc_NSD_all_eta_pt_acc_prim");
2448   aFolderObj->Add(h2D);
2449
2450   //
2451
2452   h = fMCAllPrimTrackMultHist1->Projection(0);
2453   h->SetName("mc_all_pt_acc_prim");
2454   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2455   hs->SetName("mc_all_pt_acc_prim_s");
2456   aFolderObj->Add(hs);
2457
2458   hsc = (TH1D*)hs->Clone();
2459   hsc->SetName("mc_all_pt_acc_prim_s_norm");
2460   hsc->Scale(1./fMCAllEventMultHist1->Projection(1)->Integral());
2461   aFolderObj->Add(hsc);
2462
2463   h = fMCNDEventAllPrimTrackMultHist1->Projection(0);
2464   h->SetName("mc_ND_all_pt_acc_prim");
2465   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2466   hs->SetName("mc_ND_all_pt_acc_prim_s");
2467   aFolderObj->Add(hs);
2468
2469   hsc = (TH1D*)hs->Clone();
2470   hsc->SetName("mc_ND_all_pt_acc_prim_s_norm");
2471   hsc->Scale(1./fMCAllNDEventMultHist1->Projection(1)->Integral());
2472   aFolderObj->Add(hsc);
2473
2474   h = fMCNSDEventAllPrimTrackMultHist1->Projection(0);
2475   h->SetName("mc_NSD_all_pt_acc_prim");
2476   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2477   hs->SetName("mc_NSD_all_pt_acc_prim_s");
2478   aFolderObj->Add(hs);
2479
2480   hsc = (TH1D*)hs->Clone();
2481   hsc->SetName("mc_NSD_all_pt_acc_prim_s_norm");
2482   hsc->Scale(1./fMCAllNSDEventMultHist1->Projection(1)->Integral());
2483   aFolderObj->Add(hsc);
2484
2485   h = fMCTriggerPrimTrackMultHist1->Projection(0);
2486   h->SetName("mc_trigger_all_pt_acc_prim");
2487   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2488   hs->SetName("mc_trigger_all_pt_acc_prim_s");
2489   aFolderObj->Add(hs);
2490
2491   hsc = (TH1D*)hs->Clone();
2492   hsc->SetName("mc_trigger_all_pt_acc_prim_s_norm");
2493   hsc->Scale(1./fMCTriggerMultHist1->Projection(1)->Integral());
2494   aFolderObj->Add(hsc);
2495
2496   h = fMCEventPrimTrackMultHist1->Projection(0);
2497   h->SetName("mc_all_pt_acc_trig_event_prim");
2498   hs = AlidNdPtHelper::ScaleByBinWidth(h);
2499   hs->SetName("mc_all_pt_acc_trig_event_prim_s");
2500   aFolderObj->Add(hs);
2501
2502   hsc = (TH1D*)hs->Clone();
2503   hsc->SetName("mc_all_pt_acc_trig_event_prim_s_norm");
2504   hsc->Scale(1./fMCEventMultHist1->Projection(1)->Integral());
2505   aFolderObj->Add(hsc);
2506
2507   //
2508
2509   h = fMCAllPrimTrackMultHist1->Projection(1);
2510   h->SetName("mc_all_eta_acc_prim");
2511   aFolderObj->Add(h);
2512
2513   h = fMCNDEventAllPrimTrackMultHist1->Projection(1);
2514   h->SetName("mc_ND_all_eta_acc_prim");
2515   aFolderObj->Add(h);
2516
2517   h = fMCNSDEventAllPrimTrackMultHist1->Projection(1);
2518   h->SetName("mc_NSD_all_eta_acc_prim");
2519   aFolderObj->Add(h);
2520
2521   h = fMCTriggerPrimTrackMultHist1->Projection(1);
2522   h->SetName("mc_trigger_all_eta_acc_prim");
2523   aFolderObj->Add(h);
2524
2525   h = fMCEventPrimTrackMultHist1->Projection(1);
2526   h->SetName("mc_all_eta_acc_trig_event_prim");
2527   aFolderObj->Add(h);
2528
2529   //
2530   // calculate ratios (rec / mc)
2531   //
2532   
2533   hs = (TH1*)aFolderObj->FindObject("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm");
2534   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_trig_event_track_mult_eff_cont_corrected"); 
2535   hsc->Sumw2();
2536   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_pt_acc_prim_s_norm"));
2537   aFolderObj->Add(hsc);
2538
2539   hs = (TH1*)aFolderObj->FindObject("eta_rec_trig_event_track_mult_eff_cont_corrected");
2540   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_trig_event_track_mult_eff_cont_corrected"); 
2541   hsc->Sumw2();
2542   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_eta_acc_prim"));
2543   aFolderObj->Add(hsc);
2544
2545   //
2546   hs = (TH1*)aFolderObj->FindObject("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s_norm");
2547   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_ND_trig_event_track_mult_eff_cont_corrected"); 
2548   hsc->Sumw2();
2549   hsc->Divide((TH1*)aFolderObj->FindObject("mc_ND_all_pt_acc_prim_s_norm"));
2550   aFolderObj->Add(hsc);
2551
2552   hs = (TH1*)aFolderObj->FindObject("eta_rec_ND_trig_event_track_mult_eff_cont_corrected");
2553   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_ND_trig_event_track_mult_eff_cont_corrected"); 
2554   hsc->Sumw2();
2555   hsc->Divide((TH1*)aFolderObj->FindObject("mc_ND_all_eta_acc_prim"));
2556   aFolderObj->Add(hsc);
2557
2558   //
2559   hs = (TH1*)aFolderObj->FindObject("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm");
2560   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_NSD_trig_event_track_mult_eff_cont_corrected"); 
2561   hsc->Sumw2();
2562   hsc->Divide((TH1*)aFolderObj->FindObject("mc_NSD_all_pt_acc_prim_s_norm"));
2563   aFolderObj->Add(hsc);
2564
2565   hs = (TH1*)aFolderObj->FindObject("eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");
2566   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_NSD_trig_event_track_mult_eff_cont_corrected"); 
2567   hsc->Sumw2();
2568   hsc->Divide((TH1*)aFolderObj->FindObject("mc_NSD_all_eta_acc_prim"));
2569   aFolderObj->Add(hsc);
2570
2571   //
2572   hs = (TH1*)aFolderObj->FindObject("pt_rec_event_track_mult_eff_cont_corrected_s_norm");
2573   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_event_track_mult_eff_cont_corrected"); 
2574   hsc->Sumw2();
2575   hsc->Divide((TH1*)aFolderObj->FindObject("mc_trigger_all_pt_acc_prim_s_norm"));
2576   aFolderObj->Add(hsc);
2577
2578   hs = (TH1*)aFolderObj->FindObject("eta_rec_event_track_mult_eff_cont_corrected");
2579   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_event_track_mult_eff_cont_corrected"); 
2580   hsc->Sumw2();
2581   hsc->Divide((TH1*)aFolderObj->FindObject("mc_trigger_all_eta_acc_prim"));
2582   aFolderObj->Add(hsc);
2583
2584   // track level
2585   hs = (TH1*)aFolderObj->FindObject("pt_rec_track_mult_eff_cont_corrected_s_norm");
2586   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_track_mult_eff_cont_corrected"); 
2587   hsc->Sumw2();
2588   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_pt_acc_trig_event_prim_s_norm"));
2589   aFolderObj->Add(hsc);
2590
2591   hs = (TH1*)aFolderObj->FindObject("eta_rec_track_mult_eff_cont_corrected");
2592   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_track_mult_eff_cont_corrected"); 
2593   hsc->Sumw2();
2594   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_eta_acc_trig_event_prim"));
2595   aFolderObj->Add(hsc);
2596
2597   } // end MC infor available
2598
2599   // export objects to analysis folder
2600   fCorrectionFolder = ExportToFolder(aFolderObj);
2601
2602   // delete only TObjArray
2603   if(aFolderObj) delete aFolderObj;
2604 }
2605
2606 //_____________________________________________________________________________
2607 TFolder* AlidNdPtCorrection::ExportToFolder(TObjArray * const array) 
2608 {
2609   // recreate folder avery time and export objects to new one
2610   //
2611   AlidNdPtCorrection * comp=this;
2612   TFolder *folder = comp->GetCorrectionFolder();
2613
2614   TString name, title;
2615   TFolder *newFolder = 0;
2616   Int_t i = 0;
2617   Int_t size = array->GetSize();
2618
2619   if(folder) { 
2620      // get name and title from old folder
2621      name = folder->GetName();  
2622      title = folder->GetTitle();  
2623
2624          // delete old one
2625      delete folder;
2626
2627          // create new one
2628      newFolder = CreateFolder(name.Data(),title.Data());
2629      newFolder->SetOwner();
2630
2631          // add objects to folder
2632      while(i < size) {
2633            newFolder->Add(array->At(i));
2634            i++;
2635          }
2636   }
2637
2638 return newFolder;
2639 }
2640
2641 //_____________________________________________________________________________
2642 TFolder* AlidNdPtCorrection::CreateFolder(TString name,TString title) { 
2643 // create folder for analysed histograms
2644 //
2645 TFolder *folder = 0;
2646   folder = new TFolder(name.Data(),title.Data());
2647
2648   return folder;
2649 }