removed duplications in AlidNdPtHelper and AliPWG0Helper and small changes
[u/mrichter/AliRoot.git] / PWG0 / dNdPt / AlidNdPtCorrection.cxx
1 /**************************************************************************\r
2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *\r
3  *                                                                        *\r
4  * Author: The ALICE Off-line Project.                                    *\r
5  * Contributors are mentioned in the code where appropriate.              *\r
6  *                                                                        *\r
7  * Permission to use, copy, modify and distribute this software and its   *\r
8  * documentation strictly for non-commercial purposes is hereby granted   *\r
9  * without fee, provided that the above copyright notice appears in all   *\r
10  * copies and that both the copyright notice and this permission notice   *\r
11  * appear in the supporting documentation. The authors make no claims     *\r
12  * about the suitability of this software for any purpose. It is          *\r
13  * provided "as is" without express or implied warranty.                  *\r
14  **************************************************************************/\r
15 \r
16 #include <iostream>\r
17 \r
18 #include "TFile.h"\r
19 #include "TCint.h"\r
20 #include "TH1.h"\r
21 #include "TH2.h"\r
22 \r
23 #include "AliHeader.h"  \r
24 #include "AliGenEventHeader.h"  \r
25 #include "AliStack.h"  \r
26 #include "AliESDEvent.h"  \r
27 #include "AliMCEvent.h"  \r
28 #include "AliESDtrackCuts.h"  \r
29 #include "AliLog.h" \r
30 \r
31 #include "AlidNdPtEventCuts.h"\r
32 #include "AlidNdPtAcceptanceCuts.h"\r
33 \r
34 #include "AliPWG0Helper.h"\r
35 #include "AlidNdPtHelper.h"\r
36 #include "AlidNdPtCorrection.h"\r
37 \r
38 using namespace std;\r
39 \r
40 ClassImp(AlidNdPtCorrection)\r
41 \r
42 //_____________________________________________________________________________\r
43 //AlidNdPtCorrection::AlidNdPtCorrection(): TNamed(),\r
44   AlidNdPtCorrection::AlidNdPtCorrection(): AlidNdPt(),\r
45   fCorrectionFolder(0),\r
46   fMCEventHist1(0),\r
47   fRecEventHist1(0),\r
48   fMCAllEventMultHist1(0),\r
49   fMCAllNDEventMultHist1(0),\r
50   fMCAllNSDEventMultHist1(0),\r
51   fMCTriggerMultHist1(0),\r
52   fMCEventMultHist1(0),\r
53   fMCAllPrimTrackMultHist1(0),\r
54   fMCNDEventAllPrimTrackMultHist1(0),\r
55   fMCNSDEventAllPrimTrackMultHist1(0),\r
56   fMCTriggerPrimTrackMultHist1(0),\r
57   fMCEventPrimTrackMultHist1(0),\r
58   fEventMultCorrelationMatrix(0),\r
59   fZvNorm(0),\r
60   fZvEmptyEventsNorm(0),\r
61   fCorrTriggerMBtoInelEventMatrix(0),\r
62   fCorrTriggerMBtoNDEventMatrix(0),\r
63   fCorrTriggerMBtoNSDEventMatrix(0),\r
64   fCorrEventMatrix(0),\r
65   fCorrTriggerMBtoInelTrackEventMatrix(0),\r
66   fCorrTriggerMBtoNDTrackEventMatrix(0),\r
67   fCorrTriggerMBtoNSDTrackEventMatrix(0),\r
68   fCorrTrackEventMatrix(0),\r
69   fCorrTrackMatrix(0),\r
70   fContTrackMatrix(0),\r
71   fContMultTrackMatrix(0),\r
72   fCorrMatrixFileName("")\r
73 {\r
74   // default constructor\r
75   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { \r
76     fRecTrackHist1[i]=0;     \r
77   }\r
78 \r
79   for(Int_t i=0; i<8; i++) { \r
80     fCorrRecTrackMultHist1[i] = 0;\r
81   }\r
82 \r
83   for(Int_t i=0; i<5; i++) { \r
84     fCorrRecEventHist1[i] = 0;\r
85     fCorrRecEventHist2[i] = 0;\r
86   }\r
87 \r
88   Init();\r
89 }\r
90 \r
91 //_____________________________________________________________________________\r
92 AlidNdPtCorrection::AlidNdPtCorrection(Char_t* name, Char_t* title, TString corrMatrixFileName): AlidNdPt(name,title),\r
93   fCorrectionFolder(0),\r
94   fMCEventHist1(0),\r
95   fRecEventHist1(0),\r
96   fMCAllEventMultHist1(0),\r
97   fMCAllNDEventMultHist1(0),\r
98   fMCAllNSDEventMultHist1(0),\r
99   fMCTriggerMultHist1(0),\r
100   fMCEventMultHist1(0),\r
101   fMCAllPrimTrackMultHist1(0),\r
102   fMCNDEventAllPrimTrackMultHist1(0),\r
103   fMCNSDEventAllPrimTrackMultHist1(0),\r
104   fMCTriggerPrimTrackMultHist1(0),\r
105   fMCEventPrimTrackMultHist1(0),\r
106   fEventMultCorrelationMatrix(0),\r
107   fZvNorm(0),\r
108   fZvEmptyEventsNorm(0),\r
109   fCorrTriggerMBtoInelEventMatrix(0),\r
110   fCorrTriggerMBtoNDEventMatrix(0),\r
111   fCorrTriggerMBtoNSDEventMatrix(0),\r
112   fCorrEventMatrix(0),\r
113   fCorrTriggerMBtoInelTrackEventMatrix(0),\r
114   fCorrTriggerMBtoNDTrackEventMatrix(0),\r
115   fCorrTriggerMBtoNSDTrackEventMatrix(0),\r
116   fCorrTrackEventMatrix(0),\r
117   fCorrTrackMatrix(0),\r
118   fContTrackMatrix(0),\r
119   fContMultTrackMatrix(0),\r
120   fCorrMatrixFileName(corrMatrixFileName)\r
121 {\r
122   // constructor\r
123   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { \r
124     fRecTrackHist1[i]=0;     \r
125   }\r
126 \r
127   for(Int_t i=0; i<8; i++) { \r
128     fCorrRecTrackMultHist1[i] = 0;\r
129     fPtvsPt[i] = 0;\r
130   }\r
131 \r
132   for(Int_t i=0; i<5; i++) { \r
133     fCorrRecEventHist1[i] = 0;\r
134     fCorrRecEventHist2[i] = 0;\r
135   }\r
136 \r
137   Init();\r
138 }\r
139 \r
140 //_____________________________________________________________________________\r
141 AlidNdPtCorrection::~AlidNdPtCorrection() {\r
142   // \r
143   // destructor\r
144   //\r
145   if(fMCEventHist1) delete fMCEventHist1; fMCEventHist1=0;\r
146   if(fRecEventHist1) delete fRecEventHist1; fRecEventHist1=0;\r
147 \r
148   if(fMCAllEventMultHist1) delete fMCAllEventMultHist1; fMCAllEventMultHist1=0;\r
149   if(fMCAllNDEventMultHist1) delete fMCAllNDEventMultHist1; fMCAllNDEventMultHist1=0;\r
150   if(fMCAllNSDEventMultHist1) delete fMCAllNSDEventMultHist1; fMCAllNSDEventMultHist1=0;\r
151   if(fMCTriggerMultHist1) delete fMCTriggerMultHist1; fMCTriggerMultHist1=0;\r
152   if(fMCEventMultHist1) delete fMCEventMultHist1; fMCEventMultHist1=0;\r
153 \r
154   if(fMCAllPrimTrackMultHist1) delete fMCAllPrimTrackMultHist1; fMCAllPrimTrackMultHist1=0;\r
155   if(fMCNDEventAllPrimTrackMultHist1) delete fMCNDEventAllPrimTrackMultHist1; fMCNDEventAllPrimTrackMultHist1=0;\r
156   if(fMCNSDEventAllPrimTrackMultHist1) delete fMCNSDEventAllPrimTrackMultHist1; fMCNSDEventAllPrimTrackMultHist1=0;\r
157   if(fMCTriggerPrimTrackMultHist1) delete fMCTriggerPrimTrackMultHist1; fMCTriggerPrimTrackMultHist1=0;\r
158   if(fMCEventPrimTrackMultHist1) delete fMCEventPrimTrackMultHist1; fMCEventPrimTrackMultHist1=0;\r
159 \r
160   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) { \r
161     if(fRecTrackHist1[i]) delete fRecTrackHist1[i]; fRecTrackHist1[i]=0;\r
162   }\r
163 \r
164   for(Int_t i=0; i<8; i++) { \r
165     if(fCorrRecTrackMultHist1[i]) delete fCorrRecTrackMultHist1[i]; fCorrRecTrackMultHist1[i]=0;\r
166     if(fPtvsPt[i]) delete fPtvsPt[i]; fPtvsPt[i]=0;\r
167   }\r
168 \r
169   for(Int_t i=0; i<5; i++) { \r
170     if(fCorrRecEventHist1[i]) delete fCorrRecEventHist1[i]; fCorrRecEventHist1[i]=0;\r
171     if(fCorrRecEventHist2[i]) delete fCorrRecEventHist2[i]; fCorrRecEventHist2[i]=0;\r
172   }\r
173 \r
174   if(fCorrectionFolder) delete fCorrectionFolder; fCorrectionFolder=0;\r
175 }\r
176 \r
177 //_____________________________________________________________________________\r
178 void AlidNdPtCorrection::Init(){\r
179   //\r
180   // Init histograms\r
181   //\r
182   const Int_t etaNbins = 30; \r
183   const Int_t zvNbins = 12;\r
184 \r
185   // UA1 bining\r
186   //const Int_t ptNbins = 52; \r
187   //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 }; \r
188 \r
189   const Int_t ptNbins = 56; \r
190   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};\r
191   Double_t binsEta[etaNbins+1] = {-1.5,-1.4,-1.3,-1.2,-1.1,-1.0,-0.9,-0.8,-0.7,-0.6,-0.5,-0.4,-0.3,-0.2,-0.1,0.,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,1.0,1.1,1.2,1.3,1.4,1.5};\r
192   Double_t binsZv[zvNbins+1] = {-30.,-25.,-20.,-15.,-10.,-5.,0.,5.,10.,15.,20.,25.,30.};\r
193 \r
194 \r
195   //\r
196   Int_t binsMCEventHist1[3]={100,100,140};\r
197   Double_t minMCEventHist1[3]={-0.1,-0.1,-35.}; \r
198   Double_t maxMCEventHist1[3]={0.1,0.1,35.}; \r
199   fMCEventHist1 = new THnSparseF("fMCEventHist1","mcXv:mcYv:mcZv",3,binsMCEventHist1,minMCEventHist1,maxMCEventHist1);\r
200   fMCEventHist1->GetAxis(0)->SetTitle("mcXv (cm)");\r
201   fMCEventHist1->GetAxis(1)->SetTitle("mcYv (cm)");\r
202   fMCEventHist1->GetAxis(2)->SetTitle("mcZv (cm)");\r
203   fMCEventHist1->Sumw2();\r
204 \r
205   //\r
206   Int_t binsRecEventHist1[3]={100,100,140};\r
207   Double_t minRecEventHist1[3]={-3.,-3.,-35.}; \r
208   Double_t maxRecEventHist1[3]={3.,3.,35.}; \r
209   fRecEventHist1 = new THnSparseF("fRecEventHist1","Xv:Yv:Zv",3,binsRecEventHist1,minRecEventHist1,maxRecEventHist1);\r
210   fRecEventHist1->GetAxis(0)->SetTitle("Xv (cm)");\r
211   fRecEventHist1->GetAxis(1)->SetTitle("Yv (cm)");\r
212   fRecEventHist1->GetAxis(2)->SetTitle("Zv (cm)");\r
213   fRecEventHist1->Sumw2();\r
214 \r
215   //\r
216   char name[256];\r
217   char title[256];\r
218 \r
219   Int_t binsMCAllPrimTrackMultHist1[3]={ptNbins,etaNbins,150};\r
220   Double_t minMCAllPrimTrackMultHist1[3]={0.,-1.,-0.5}; \r
221   Double_t maxMCAllPrimTrackMultHist1[3]={20.,1.,149.5}; \r
222   sprintf(name,"fMCAllPrimTrackMultHist1");\r
223   sprintf(title,"mcPt:mcEta:multiplicity");\r
224   \r
225   fMCAllPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCAllPrimTrackMultHist1,minMCAllPrimTrackMultHist1,maxMCAllPrimTrackMultHist1);\r
226   fMCAllPrimTrackMultHist1->SetBinEdges(0,binsPt);\r
227   fMCAllPrimTrackMultHist1->SetBinEdges(1,binsEta);\r
228   fMCAllPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
229   fMCAllPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");\r
230   fMCAllPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");\r
231   fMCAllPrimTrackMultHist1->Sumw2();\r
232 \r
233   Int_t binsMCNDEventAllPrimTrackMultHist1[3]={ptNbins,etaNbins,150};\r
234   Double_t minMCNDEventAllPrimTrackMultHist1[3]={0.,-1.,-0.5}; \r
235   Double_t maxMCNDEventAllPrimTrackMultHist1[3]={20.,1.,149.5}; \r
236   sprintf(name,"fMCNDEventAllPrimTrackMultHist1");\r
237   sprintf(title,"mcPt:mcEta:multiplicity");\r
238   \r
239   fMCNDEventAllPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCNDEventAllPrimTrackMultHist1,minMCNDEventAllPrimTrackMultHist1,maxMCNDEventAllPrimTrackMultHist1);\r
240   fMCNDEventAllPrimTrackMultHist1->SetBinEdges(0,binsPt);\r
241   fMCNDEventAllPrimTrackMultHist1->SetBinEdges(1,binsEta);\r
242   fMCNDEventAllPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
243   fMCNDEventAllPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");\r
244   fMCNDEventAllPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");\r
245   fMCNDEventAllPrimTrackMultHist1->Sumw2();\r
246 \r
247   Int_t binsMCNSDEventAllPrimTrackMultHist1[3]={ptNbins,etaNbins,150};\r
248   Double_t minMCNSDEventAllPrimTrackMultHist1[3]={0.,-1.,-0.5}; \r
249   Double_t maxMCNSDEventAllPrimTrackMultHist1[3]={20.,1.,140.9}; \r
250   sprintf(name,"fMCNSDEventAllPrimTrackMultHist1");\r
251   sprintf(title,"mcPt:mcEta:multiplicity");\r
252   \r
253   fMCNSDEventAllPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCNSDEventAllPrimTrackMultHist1,minMCNSDEventAllPrimTrackMultHist1,maxMCNSDEventAllPrimTrackMultHist1);\r
254   fMCNSDEventAllPrimTrackMultHist1->SetBinEdges(0,binsPt);\r
255   fMCNSDEventAllPrimTrackMultHist1->SetBinEdges(1,binsEta);\r
256   fMCNSDEventAllPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
257   fMCNSDEventAllPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");\r
258   fMCNSDEventAllPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");\r
259   fMCNSDEventAllPrimTrackMultHist1->Sumw2();\r
260 \r
261   Int_t binsMCEventTriggerPrimTrackMultHist1[3]={ptNbins,etaNbins,150};\r
262   Double_t minMCEventTriggerPrimTrackMultHist1[3]={0.,-1.,-0.5}; \r
263   Double_t maxMCEventTriggerPrimTrackMultHist1[3]={20.,1.,149.5}; \r
264   sprintf(name,"fMCTriggerPrimTrackMultHist1");\r
265   sprintf(title,"mcPt:mcEta:multiplicity");\r
266   \r
267   fMCTriggerPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCEventTriggerPrimTrackMultHist1,minMCEventTriggerPrimTrackMultHist1,maxMCEventTriggerPrimTrackMultHist1);\r
268   fMCTriggerPrimTrackMultHist1->SetBinEdges(0,binsPt);\r
269   fMCTriggerPrimTrackMultHist1->SetBinEdges(1,binsEta);\r
270   fMCTriggerPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
271   fMCTriggerPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");\r
272   fMCTriggerPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");\r
273   fMCTriggerPrimTrackMultHist1->Sumw2();\r
274 \r
275   Int_t binsMCEventPrimTrackMultHist1[3]={ptNbins,etaNbins,150};\r
276   Double_t minMCEventPrimTrackMultHist1[3]={0.,-1.,-0.5}; \r
277   Double_t maxMCEventPrimTrackMultHist1[3]={20.,1.,149.5}; \r
278   sprintf(name,"fMCEventPrimTrackMultHist1");\r
279   sprintf(title,"mcPt:mcEta:multiplicity");\r
280   \r
281   fMCEventPrimTrackMultHist1 = new THnSparseF(name,title,3,binsMCEventPrimTrackMultHist1,minMCEventPrimTrackMultHist1,maxMCEventPrimTrackMultHist1);\r
282   fMCEventPrimTrackMultHist1->SetBinEdges(0,binsPt);\r
283   fMCEventPrimTrackMultHist1->SetBinEdges(1,binsEta);\r
284   fMCEventPrimTrackMultHist1->GetAxis(0)->SetTitle("mcPt (GeV/c)");\r
285   fMCEventPrimTrackMultHist1->GetAxis(1)->SetTitle("mcEta");\r
286   fMCEventPrimTrackMultHist1->GetAxis(2)->SetTitle("multiplicity");\r
287   fMCEventPrimTrackMultHist1->Sumw2();\r
288 \r
289   for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) \r
290   {\r
291     // THnSparse track histograms\r
292     //\r
293     Int_t binsRecTrackHist1[3]={ptNbins,etaNbins,90};\r
294     Double_t minRecTrackHist1[3]={0.,-1.,0.}; \r
295     Double_t maxRecTrackHist1[3]={10.,1.,2.*TMath::Pi()};\r
296     sprintf(name,"fRecTrackHist1_%d",i);\r
297     sprintf(title,"Pt:Eta:Phi");\r
298   \r
299     fRecTrackHist1[i] = new THnSparseF(name,title,3,binsRecTrackHist1,minRecTrackHist1,maxRecTrackHist1);\r
300     fRecTrackHist1[i]->SetBinEdges(0,binsPt);\r
301     fRecTrackHist1[i]->SetBinEdges(1,binsEta);\r
302     fRecTrackHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");\r
303     fRecTrackHist1[i]->GetAxis(1)->SetTitle("Eta");\r
304     fRecTrackHist1[i]->GetAxis(2)->SetTitle("Phi (rad)");\r
305     fRecTrackHist1[i]->Sumw2();\r
306   }\r
307 \r
308   Int_t binsCorrRecTrackMultHist1[3]={ptNbins,etaNbins,150};\r
309   Double_t minCorrRecTrackMultHist1[3]={0.,-1.,-0.5}; \r
310   Double_t maxCorrRecTrackMultHist1[3]={20.,1.,149.5};\r
311 \r
312   Int_t binsPtvsPt[3]={ptNbins,320};\r
313   Double_t minPtvsPt[3]={0.,0.}; \r
314   Double_t maxPtvsPt[3]={20.,16.};\r
315 \r
316   for(Int_t i=0; i<8; i++) \r
317   {\r
318     // THnSparse track histograms\r
319     sprintf(name,"fCorrRecTrackMultHist1_%d",i);\r
320     sprintf(title,"Pt:Eta:mult");\r
321     fCorrRecTrackMultHist1[i] = new THnSparseF(name,title,3,binsCorrRecTrackMultHist1,minCorrRecTrackMultHist1,maxCorrRecTrackMultHist1);\r
322     fCorrRecTrackMultHist1[i]->SetBinEdges(0,binsPt);\r
323     fCorrRecTrackMultHist1[i]->SetBinEdges(1,binsEta);\r
324     fCorrRecTrackMultHist1[i]->GetAxis(0)->SetTitle("Pt (GeV/c)");\r
325     fCorrRecTrackMultHist1[i]->GetAxis(1)->SetTitle("Eta");\r
326     fCorrRecTrackMultHist1[i]->GetAxis(2)->SetTitle("multiplicity");\r
327     fCorrRecTrackMultHist1[i]->Sumw2();\r
328 \r
329     sprintf(name,"fPtvsPt_%d",i);\r
330     sprintf(title,"Pt:Pt");\r
331     fPtvsPt[i] = new THnSparseF(name,title,2,binsPtvsPt,minPtvsPt,maxPtvsPt);\r
332     fPtvsPt[i]->SetBinEdges(0,binsPt);\r
333     fPtvsPt[i]->GetAxis(0)->SetTitle("Pt (GeV/c)"); \r
334     fPtvsPt[i]->GetAxis(1)->SetTitle("Pt (GeV/c)"); \r
335     fPtvsPt[i]->Sumw2();\r
336   }\r
337 \r
338   Int_t binsEventMatrix[2]={zvNbins,150};\r
339   Double_t minEventMatrix[2]={-25.,-0.5};\r
340   Double_t maxEventMatrix[2]={25.,149.5};\r
341 \r
342   fMCAllEventMultHist1 = new THnSparseF("fMCAllEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
343   fMCAllEventMultHist1->SetBinEdges(0,binsZv);\r
344   fMCAllEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");\r
345   fMCAllEventMultHist1->GetAxis(1)->SetTitle("multiplicity");\r
346   fMCAllEventMultHist1->Sumw2();\r
347 \r
348   fMCAllNDEventMultHist1 = new THnSparseF("fMCAllNDEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
349   fMCAllNDEventMultHist1->SetBinEdges(0,binsZv);\r
350   fMCAllNDEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");\r
351   fMCAllNDEventMultHist1->GetAxis(1)->SetTitle("multiplicity");\r
352   fMCAllNDEventMultHist1->Sumw2();\r
353 \r
354   fMCAllNSDEventMultHist1 = new THnSparseF("fMCAllNSDEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
355   fMCAllNSDEventMultHist1->SetBinEdges(0,binsZv);\r
356   fMCAllNSDEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");\r
357   fMCAllNSDEventMultHist1->GetAxis(1)->SetTitle("multiplicity");\r
358   fMCAllNSDEventMultHist1->Sumw2();\r
359 \r
360   fMCTriggerMultHist1 = new THnSparseF("fMCTriggerMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
361   fMCTriggerMultHist1->SetBinEdges(0,binsZv);\r
362   fMCTriggerMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");\r
363   fMCTriggerMultHist1->GetAxis(1)->SetTitle("multiplicity");\r
364   fMCTriggerMultHist1->Sumw2();\r
365 \r
366   fMCEventMultHist1 = new THnSparseF("fMCEventMultHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
367   fMCEventMultHist1->SetBinEdges(0,binsZv);\r
368   fMCEventMultHist1->GetAxis(0)->SetTitle("mcZv (cm)");\r
369   fMCEventMultHist1->GetAxis(1)->SetTitle("multiplicity");\r
370   fMCEventMultHist1->Sumw2();\r
371 \r
372   for(Int_t i=0; i<5; i++) \r
373   {\r
374     // event corrected histograms\r
375     sprintf(name,"fCorrRecEventHist1_%d",i);\r
376     sprintf(title,"mcZv:mult");\r
377     fCorrRecEventHist1[i] = new THnSparseF("fCorrRecEventHist1","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
378     fCorrRecEventHist1[i]->SetBinEdges(0,binsZv);\r
379     fCorrRecEventHist1[i]->GetAxis(0)->SetTitle("Zv (cm)");\r
380     fCorrRecEventHist1[i]->GetAxis(1)->SetTitle("multiplicity");\r
381     fCorrRecEventHist1[i]->Sumw2();\r
382 \r
383     // empty event corrected histograms\r
384     sprintf(name,"fCorrRecEventHist2_%d",i);\r
385     sprintf(title,"mcZv:mult");\r
386     fCorrRecEventHist2[i] = new THnSparseF("fCorrRecEventHist2","mcZv:mult",2,binsEventMatrix,minEventMatrix,maxEventMatrix);\r
387     fCorrRecEventHist2[i]->SetBinEdges(0,binsZv);\r
388     fCorrRecEventHist2[i]->GetAxis(0)->SetTitle("Zv (cm)");\r
389     fCorrRecEventHist2[i]->GetAxis(1)->SetTitle("multiplicity");\r
390     fCorrRecEventHist2[i]->Sumw2();\r
391   }\r
392 \r
393   // init output folder\r
394   fCorrectionFolder = CreateFolder("folderdNdPt","Correction dNdPt Folder");\r
395 \r
396   // init correction matrices\r
397   TFile *file = (TFile *)TFile::Open(fCorrMatrixFileName.Data()); \r
398   if(!file) { \r
399     AliDebug(AliLog::kError, "file with efficiency matrices not available");\r
400     printf("file with efficiency matrices not available \n");\r
401   } else {\r
402     TFolder *folder = (TFolder *)file->FindObjectAny("folderdNdPt");\r
403     if(!folder) { \r
404       AliDebug(AliLog::kError, "file without folderdNdPt");\r
405       printf("file without folderdNdPt \n");\r
406     } else {\r
407       // rec. event mult vs true multiplicity \r
408       fEventMultCorrelationMatrix = (THnSparseF*)folder->FindObject("event_mult_correlation_matrix");\r
409       if(!fEventMultCorrelationMatrix) {\r
410          Printf("No %s matrix \n", "event_mult_correlation_matrix");\r
411       }\r
412 \r
413       //\r
414       // event level corrections (zv,mult_MB)\r
415       //\r
416  \r
417       // trigger bias correction (MBtoND) \r
418       fCorrTriggerMBtoNDEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_trig_MBtoND_corr_matrix");\r
419       if(!fCorrTriggerMBtoNDEventMatrix) {\r
420          Printf("No %s matrix \n", "zv_mult_trig_MBtoND_corr_matrix");\r
421       }\r
422 \r
423       // trigger bias correction (MBtoNSD)\r
424       fCorrTriggerMBtoNSDEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_trig_MBtoNSD_corr_matrix");\r
425       if(!fCorrTriggerMBtoNSDEventMatrix) {\r
426          Printf("No %s matrix \n", "zv_mult_trig_MBtoNSD_corr_matrix");\r
427       }\r
428 \r
429       // trigger bias correction (MBtoInel)\r
430       fCorrTriggerMBtoInelEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_trig_MBtoInel_corr_matrix");\r
431       if(!fCorrTriggerMBtoInelEventMatrix) {\r
432          Printf("No %s matrix \n", "zv_mult_trig_MBtoInel_corr_matrix"); \r
433       }\r
434      \r
435       // vertex reconstruction efficiency correction\r
436       fCorrEventMatrix = (THnSparseF*)folder->FindObject("zv_mult_event_corr_matrix");\r
437       if(!fCorrEventMatrix) {\r
438          Printf("No %s matrix \n", "zv_mult_event_corr_matrix");\r
439       }\r
440 \r
441       //\r
442       // histogram needed for empty events corrections\r
443       //\r
444       fZvNorm = (TH1D*)folder->FindObject("zv_distribution_norm");\r
445       if(!fZvNorm) {\r
446          Printf("No %s matrix \n", "fZvNorm");\r
447       }\r
448 \r
449       fZvEmptyEventsNorm = (TH1D*)folder->FindObject("zv_empty_events_norm");\r
450       if(!fZvEmptyEventsNorm) {\r
451          Printf("No %s matrix \n", "fZvEmptyEventsNorm");\r
452       }\r
453 \r
454       //\r
455       // track-event level corrections (zv,pt,eta)\r
456       //\r
457 \r
458       // trigger bias correction (MBtoND) \r
459       fCorrTriggerMBtoNDTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_trig_MBtoND_corr_matrix");\r
460       if(!fCorrTriggerMBtoNDTrackEventMatrix) {\r
461          Printf("No %s matrix \n", "zv_pt_eta_track_trig_MBtoND_corr_matrix");\r
462       }\r
463 \r
464       // trigger bias correction (MBtoNSD)\r
465       fCorrTriggerMBtoNSDTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_trig_MBtoNSD_corr_matrix");\r
466       if(!fCorrTriggerMBtoNSDTrackEventMatrix) {\r
467          Printf("No %s matrix \n", "zv_pt_eta_track_trig_MBtoNSD_corr_matrix");\r
468       }\r
469 \r
470       // trigger bias correction (MBtoInel) \r
471       fCorrTriggerMBtoInelTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_trig_MBtoInel_corr_matrix");\r
472       if(!fCorrTriggerMBtoInelTrackEventMatrix) {\r
473          Printf("No %s matrix \n", "zv_pt_eta_track_trig_MBtoInel_corr_matrix");\r
474       }\r
475     \r
476       // vertex reconstruction efficiency correction (zv,pt,eta)\r
477       fCorrTrackEventMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_event_corr_matrix");\r
478       if(!fCorrTrackEventMatrix) {\r
479          Printf("No %s matrix \n", "zv_pt_eta_track_event_corr_matrix");\r
480       }\r
481 \r
482       // track reconstruction efficiency correction (zv,pt,eta)\r
483       fCorrTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_corr_matrix");\r
484       if(!fCorrTrackMatrix) {\r
485          Printf("No %s matrix \n", "zv_pt_eta_track_corr_matrix");\r
486       }\r
487 \r
488       // secondary tracks contamination correction (zv,pt,eta)\r
489       fContTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_track_cont_matrix");\r
490       if(!fContTrackMatrix) {\r
491          Printf("No %s matrix \n", "zv_pt_eta_track_cont_matrix");\r
492       }\r
493 \r
494       // multiply reconstructed tracks correction\r
495       fContMultTrackMatrix = (THnSparseF*)folder->FindObject("zv_pt_eta_mult_track_cont_matrix");\r
496       if(!fContMultTrackMatrix) {\r
497          Printf("No %s matrix \n", "zv_pt_eta_mult_track_cont_matrix");\r
498       }\r
499     }\r
500   }\r
501 \r
502 }\r
503 \r
504 //_____________________________________________________________________________\r
505 void AlidNdPtCorrection::Process(AliESDEvent *esdEvent, AliMCEvent *mcEvent)\r
506 {\r
507   //\r
508   // Process real and/or simulated events\r
509   //\r
510   if(!esdEvent) {\r
511     AliDebug(AliLog::kError, "esdEvent not available");\r
512     return;\r
513   }\r
514 \r
515   // get selection cuts\r
516   AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
517   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
518   AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
519 \r
520   if(!evtCuts || !accCuts  || !esdTrackCuts) {\r
521     AliDebug(AliLog::kError, "cuts not available");\r
522     return;\r
523   }\r
524 \r
525   // trigger selection\r
526   Bool_t isEventTriggered = kTRUE;\r
527   if(evtCuts->IsTriggerRequired())  {\r
528     isEventTriggered = AliPWG0Helper::IsEventTriggered(esdEvent->GetTriggerMask(), GetTrigger());\r
529   }\r
530 \r
531   // use MC information\r
532   AliHeader* header = 0;\r
533   AliGenEventHeader* genHeader = 0;\r
534   AliStack* stack = 0;\r
535   TArrayF vtxMC(3);\r
536   AliPWG0Helper::MCProcessType evtType = AliPWG0Helper::kInvalidProcess;\r
537   Int_t multMCTrueTracks = 0;\r
538 \r
539   if(IsUseMCInfo())\r
540   {\r
541     if(!mcEvent) {\r
542       AliDebug(AliLog::kError, "mcEvent not available");\r
543       return;\r
544     }\r
545 \r
546     // get MC event header\r
547     header = mcEvent->Header();\r
548     if (!header) {\r
549       AliDebug(AliLog::kError, "Header not available");\r
550       return;\r
551     }\r
552     // MC particle stack\r
553     stack = mcEvent->Stack();\r
554     if (!stack) {\r
555       AliDebug(AliLog::kError, "Stack not available");\r
556       return;\r
557     }\r
558 \r
559     // get event type (ND=0x1, DD=0x2, SD=0x4)\r
560     evtType = AliPWG0Helper::GetEventProcessType(header);\r
561     //Printf("evtType %d \n", evtType);\r
562     AliDebug(AliLog::kDebug+1, Form("Found process type %d", evtType));\r
563 \r
564     // get MC vertex\r
565     genHeader = header->GenEventHeader();\r
566     if (!genHeader) {\r
567       AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");\r
568       return;\r
569     }\r
570     genHeader->PrimaryVertex(vtxMC);\r
571 \r
572     // Fill MC event histogram\r
573     Double_t vMCEventHist1[3]={vtxMC[0],vtxMC[1],vtxMC[2]};\r
574     fMCEventHist1->Fill(vMCEventHist1);\r
575 \r
576     // multipliticy of all MC primary tracks\r
577     // in Zvtx, pt and eta ranges\r
578     multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);\r
579 \r
580   } // end bUseMC\r
581 \r
582   // get reconstructed vertex  \r
583   const AliESDVertex* vtxESD = 0; \r
584   Bool_t isRecVertex = kFALSE;\r
585   if(evtCuts->IsRecVertexRequired()) \r
586   {\r
587     Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();\r
588     Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();\r
589     vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints); \r
590     isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, GetAnalysisMode(), kFALSE); // should be moved to AcceptEvent\r
591   }\r
592   if( IsUseMCInfo() &&  !evtCuts->IsRecVertexRequired() ) {\r
593     vtxESD = new AliESDVertex(vtxMC[2],10.,genHeader->NProduced());\r
594     isRecVertex = kTRUE;\r
595   }\r
596   Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex; \r
597   //printf("isEventOK %d \n",isEventOK);\r
598 \r
599   //\r
600   // get multiplicity of min. bias tracks\r
601   //\r
602   Int_t multMBTracks = 0; \r
603   if(GetAnalysisMode() == AlidNdPtHelper::kTPC || GetAnalysisMode() == AlidNdPtHelper::kMCRec) {  \r
604     multMBTracks = AlidNdPtHelper::GetTPCMBTrackMult(esdEvent,evtCuts,accCuts,esdTrackCuts);\r
605   } \r
606   else if(GetAnalysisMode() == AlidNdPtHelper::kTPCSPDvtx || GetAnalysisMode() == AlidNdPtHelper::kMCRec) {\r
607     multMBTracks = AlidNdPtHelper::GetSPDMBTrackMult(esdEvent,0.0);\r
608   } \r
609   else {\r
610     AliDebug(AliLog::kError, Form("Found analysis type %s", GetAnalysisMode()));\r
611     return; \r
612   }\r
613 \r
614   //\r
615   // correct event and track histograms\r
616   //\r
617   TObjArray *allChargedTracks=0;\r
618   Int_t multRec=0, multRecTemp=0;\r
619   Int_t *labelsRec=0;\r
620 \r
621   if(isEventOK && isEventTriggered)\r
622   {\r
623     // get all charged tracks\r
624     //allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,vtxESD,GetAnalysisMode());\r
625     allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());\r
626     if(!allChargedTracks) return;\r
627 \r
628     Int_t entries = allChargedTracks->GetEntries();\r
629     labelsRec = new Int_t[entries];\r
630 \r
631     // calculate mult of reconstructed tracks\r
632     for(Int_t i=0; i<entries;++i) \r
633     {\r
634       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
635       if(!track) continue;\r
636       if(accCuts->AcceptTrack(track)) \r
637       {\r
638         if(esdTrackCuts->AcceptTrack(track)) \r
639         {\r
640           multRecTemp++;\r
641         }\r
642       }\r
643     }  \r
644 \r
645     for(Int_t i=0; i<entries;++i) \r
646     {\r
647       AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
648       if(!track) continue;\r
649         \r
650       if(accCuts->AcceptTrack(track)) \r
651       {\r
652         if(esdTrackCuts->AcceptTrack(track)) \r
653         {\r
654           // track-level corrections\r
655           if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) { \r
656             FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxMC[2],multRecTemp); \r
657           } else {\r
658             FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxESD->GetZv(),multRecTemp); \r
659           }\r
660           labelsRec[multRec] = TMath::Abs(track->GetLabel());\r
661           multRec++;\r
662         }\r
663       }\r
664     }\r
665     // event-level corrections\r
666     if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) { \r
667       FillHistograms(AlidNdPtHelper::kRecEvents,vtxMC[2],multMBTracks);\r
668     }\r
669     else {\r
670       FillHistograms(AlidNdPtHelper::kRecEvents,vtxESD->GetZv(),multMBTracks);\r
671     }\r
672 \r
673     // control event histograms\r
674     Double_t vRecEventHist1[3] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv()};\r
675     fRecEventHist1->Fill(vRecEventHist1);\r
676   } \r
677 \r
678   // empty events corrections\r
679   // no reconstructed zv\r
680   if(isEventTriggered && multMBTracks==0) {\r
681     if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) {\r
682       FillHistograms(AlidNdPtHelper::kTriggeredEvents,vtxMC[2],multMBTracks);\r
683     }\r
684     else {\r
685       Double_t zv = fZvNorm->GetRandom();\r
686       FillHistograms(AlidNdPtHelper::kTriggeredEvents,zv,multMBTracks);\r
687     }\r
688   }\r
689 \r
690   if(IsUseMCInfo())  \r
691   {\r
692     // select MC events \r
693     if(evtCuts->AcceptMCEvent(mcEvent))\r
694     {\r
695       //\r
696       // event histograms\r
697       //\r
698       Double_t vMCEventMatrix[2] = {vtxMC[2],multMBTracks};\r
699       fMCAllEventMultHist1->Fill(vMCEventMatrix);\r
700       if(evtType == AliPWG0Helper::kND) {\r
701         fMCAllNDEventMultHist1->Fill(vMCEventMatrix);\r
702       }\r
703       if(evtType != AliPWG0Helper::kSD) {\r
704         fMCAllNSDEventMultHist1->Fill(vMCEventMatrix);\r
705       }\r
706       if(isEventTriggered) fMCTriggerMultHist1->Fill(vMCEventMatrix);\r
707       if(isEventTriggered && isEventOK) fMCEventMultHist1->Fill(vMCEventMatrix);\r
708 \r
709       //\r
710       // MC histograms for efficiency studies\r
711       //\r
712       Int_t nPart  = stack->GetNtrack();\r
713       for (Int_t iMc = 0; iMc < nPart; ++iMc) \r
714       {\r
715         // print MC stack info\r
716         //AlidNdPtHelper::PrintMCInfo(stack,iMc);\r
717 \r
718         TParticle* particle = stack->Particle(iMc);\r
719         if (!particle)\r
720         continue;\r
721 \r
722         // only charged particles\r
723         Double_t charge = particle->GetPDG()->Charge()/3.;\r
724         if (charge == 0.0)\r
725           continue;\r
726       \r
727         // physical primary\r
728         Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
729 \r
730         // all primaries in acceptance\r
731         if(!accCuts->AcceptTrack(particle)) continue;\r
732         if(!prim) continue;\r
733 \r
734         Double_t gpt = particle->Pt();\r
735         Double_t geta = particle->Eta();\r
736 \r
737         Double_t valMCAllTrackMultHist1[3] = {gpt,geta,multRec};          \r
738         fMCAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
739         if(evtType == AliPWG0Helper::kND) {\r
740           fMCNDEventAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
741         }\r
742         if(evtType != AliPWG0Helper::kSD) {\r
743           fMCNSDEventAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
744         }\r
745         if(isEventTriggered) fMCTriggerPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
746         if(isEventTriggered && isEventOK) fMCEventPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
747       }\r
748     }\r
749   } // end bUseMC\r
750 \r
751   if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;\r
752   if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
753 \r
754   if(!evtCuts->IsRecVertexRequired() && vtxESD != NULL) \r
755   delete vtxESD;\r
756 }\r
757 \r
758 //_____________________________________________________________________________\r
759 void AlidNdPtCorrection::FillHistograms(AlidNdPtHelper::EventObject eventObj, Double_t zv, Int_t multMBTracks)\r
760 {\r
761   //\r
762   // Fill corrected histograms\r
763   //\r
764 \r
765   Double_t vEventMatrix[2] = {zv,multMBTracks};\r
766   //\r
767   // Correct for efficiency \r
768   //\r
769   if(eventObj == AlidNdPtHelper::kRecEvents && multMBTracks>0)  \r
770   {\r
771     Double_t corrToMBF = GetCorrFactZvMult(fCorrEventMatrix,zv,multMBTracks);\r
772     Double_t corrToInelF = GetCorrFactZvMult(fCorrTriggerMBtoInelEventMatrix,zv,multMBTracks);\r
773     Double_t corrToNDF = GetCorrFactZvMult(fCorrTriggerMBtoNDEventMatrix,zv,multMBTracks);\r
774     Double_t corrToNSDF = GetCorrFactZvMult(fCorrTriggerMBtoNSDEventMatrix,zv,multMBTracks);\r
775     //printf("corrToMBF %f, corrToInelF %f, corrToNDF %f corrToNSDF %f \n",corrToMBF,corrToInelF,corrToNDF,corrToNSDF);\r
776 \r
777     fCorrRecEventHist1[0]->Fill(vEventMatrix);\r
778     fCorrRecEventHist1[1]->Fill(vEventMatrix,corrToMBF);\r
779     fCorrRecEventHist1[2]->Fill(vEventMatrix,corrToMBF*corrToInelF);\r
780     fCorrRecEventHist1[3]->Fill(vEventMatrix,corrToMBF*corrToNDF);\r
781     fCorrRecEventHist1[4]->Fill(vEventMatrix,corrToMBF*corrToNSDF);\r
782   }\r
783 \r
784   if(eventObj==AlidNdPtHelper::kTriggeredEvents && multMBTracks==0) // empty triggered events\r
785   {\r
786     Int_t bin = fZvEmptyEventsNorm->FindBin(zv);\r
787     Double_t Fz = fZvEmptyEventsNorm->GetBinContent(bin);\r
788     Double_t corrToInelF0 = GetCorrFactZvMult(fCorrTriggerMBtoInelEventMatrix,zv,multMBTracks);\r
789     Double_t corrToNDF0 = GetCorrFactZvMult(fCorrTriggerMBtoNDEventMatrix,zv,multMBTracks);\r
790     Double_t corrToNSDF0 = GetCorrFactZvMult(fCorrTriggerMBtoNSDEventMatrix,zv,multMBTracks);\r
791     //printf("Fz %f, corrToInelF0 %f, corrToNDF0 %f, corrToNSDF0 %f \n",Fz,corrToInelF0,corrToNDF0,corrToNSDF0);\r
792 \r
793     fCorrRecEventHist2[0]->Fill(vEventMatrix);\r
794     fCorrRecEventHist2[1]->Fill(vEventMatrix,Fz);\r
795     fCorrRecEventHist2[2]->Fill(vEventMatrix,Fz*corrToInelF0);\r
796     fCorrRecEventHist2[3]->Fill(vEventMatrix,Fz*corrToNDF0);\r
797     fCorrRecEventHist2[4]->Fill(vEventMatrix,Fz*corrToNSDF0);\r
798   }\r
799 }\r
800 \r
801 //_____________________________________________________________________________\r
802 void AlidNdPtCorrection::FillHistograms(AliESDtrack *esdTrack, AliStack *stack, AlidNdPtHelper::TrackObject trackObj, Double_t zv, Int_t mult)\r
803 {\r
804   //\r
805   // Fill ESD track and MC histograms \r
806   //\r
807   if(!esdTrack) return;\r
808 \r
809   //Float_t q = esdTrack->Charge();\r
810   Float_t pt = esdTrack->Pt();\r
811   Float_t eta = esdTrack->Eta();\r
812   Float_t phi = esdTrack->Phi();\r
813 \r
814   if(stack && GetAnalysisMode() == AlidNdPtHelper::kMCRec) \r
815   {\r
816     Int_t label = TMath::Abs(esdTrack->GetLabel());\r
817    \r
818     TParticle* particle = stack->Particle(label);\r
819     if(!particle) return;\r
820    \r
821     Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3\r
822     if(gq==0) return;\r
823     Float_t gpt = particle->Pt();\r
824     Float_t geta = particle->Eta();\r
825     Float_t gphi = particle->Phi();\r
826 \r
827     // replace reconstructed values with MC\r
828     pt = gpt;\r
829     eta = geta;\r
830     phi = gphi;\r
831   }\r
832 \r
833   //\r
834   // Fill histograms\r
835   //\r
836   Double_t values[3] = {pt,eta,phi};      \r
837   fRecTrackHist1[trackObj]->Fill(values);\r
838 \r
839   //\r
840   // Correct for contamination and efficiency \r
841   //\r
842   if(trackObj == AlidNdPtHelper::kRecTracks || GetAnalysisMode() == AlidNdPtHelper::kMCRec)  \r
843   {\r
844     // track level corrections\r
845     Double_t trackEffF = GetCorrFactZvPtEta(fCorrTrackMatrix,zv,pt,eta);\r
846     Double_t trackContF = GetContFactZvPtEta(fContTrackMatrix,zv,pt,eta);\r
847     Double_t multTrackContF = GetContFactZvPtEta(fContMultTrackMatrix,zv,pt,eta);\r
848     //printf("zv %f, pt %f, eta %f \n",zv,pt,eta);\r
849     //printf("trackEffF %f, trackContF %f, multTrackContF %f \n", trackEffF, trackContF, multTrackContF);\r
850    \r
851     // track-event level corrections\r
852     Double_t vertexEffF = GetCorrFactZvPtEta(fCorrTrackEventMatrix,zv,pt,eta);\r
853     Double_t trigMBToInel = GetCorrFactZvPtEta(fCorrTriggerMBtoInelTrackEventMatrix,zv,pt,eta);  \r
854     Double_t trigMBToND = GetCorrFactZvPtEta(fCorrTriggerMBtoNDTrackEventMatrix,zv,pt,eta);\r
855     Double_t trigMBToNSD = GetCorrFactZvPtEta(fCorrTriggerMBtoNSDTrackEventMatrix,zv,pt,eta);\r
856     //printf("vertexEffF %f, trigMBToInel %f, trigMBToNSD %f \n", vertexEffF, trigMBToInel, trigMBToNSD);\r
857     \r
858     Double_t corrF[8] = { 1.0, \r
859                           trackContF,\r
860                           trackContF*trackEffF,\r
861                           trackContF*trackEffF*multTrackContF,\r
862                           trackContF*trackEffF*multTrackContF*vertexEffF,\r
863                           trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToInel,\r
864                           trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToND,\r
865                           trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToNSD\r
866                          }; \r
867  \r
868     // Fill histograms\r
869     Double_t valCorrRecTrackMultHist1[3] = {pt,eta,mult};         \r
870     Double_t valPtvsPt[2] = {pt,pt};      \r
871     for(Int_t i=0; i<8; i++) {\r
872       fCorrRecTrackMultHist1[i]->Fill(valCorrRecTrackMultHist1,corrF[i]);\r
873       fPtvsPt[i]->Fill(valPtvsPt,corrF[i]);\r
874     }\r
875   }\r
876 }\r
877 \r
878 //_____________________________________________________________________________\r
879 void AlidNdPtCorrection::FillHistograms(AliStack *stack, Int_t /*label*/, AlidNdPtHelper::TrackObject /*trackObj*/, Int_t /*mult*/)\r
880 {\r
881   // Fill MC histograms\r
882   if(!stack) return;\r
883 \r
884   /*\r
885   TParticle* particle = stack->Particle(label);\r
886   if(!particle) return;\r
887 \r
888   Int_t mother_pdg = -1;\r
889   TParticle* mother = 0;\r
890 \r
891   //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);\r
892   Int_t motherLabel = particle->GetMother(0); \r
893   if(motherLabel>0) mother = stack->Particle(motherLabel);\r
894   if(mother) mother_pdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only\r
895   Int_t mech = particle->GetUniqueID(); // production mechanism\r
896 \r
897   Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 \r
898   Float_t gpt = particle->Pt();\r
899   Float_t qgpt = particle->Pt() * gq;\r
900   Float_t geta = particle->Eta();\r
901   Float_t gphi = particle->Phi();\r
902   Float_t gpz = particle->Pz();\r
903 \r
904   Bool_t prim = stack->IsPhysicalPrimary(label);\r
905   Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();\r
906 \r
907   Int_t pid=-1;\r
908   if (TMath::Abs(particle->GetPdgCode()) == kElectron)         { pid = 0; }\r
909     else if (TMath::Abs(particle->GetPdgCode()) == kMuonMinus) { pid = 1; }\r
910     else if (TMath::Abs(particle->GetPdgCode()) == kPiPlus)    { pid = 2; }\r
911     else if (TMath::Abs(particle->GetPdgCode()) == kKPlus)     { pid = 3; }\r
912     else if (TMath::Abs(particle->GetPdgCode()) == kProton)    { pid = 4; }\r
913     else                                                       { pid = 5; }\r
914     */\r
915 \r
916   //if(!prim) printf("prim_mother %d, mother %d, particle %d, production mech %d\n",prim_mother->GetPdgCode(),mother->GetPdgCode(), particle->GetPdgCode(),mech);\r
917   \r
918 }\r
919 \r
920 //_____________________________________________________________________________\r
921 Double_t AlidNdPtCorrection::GetCorrFactZvPtEta(THnSparse *hist, Double_t zv, Double_t pt, Double_t eta) const {\r
922 // return correction factor F(zv,pt,eta)\r
923 \r
924  if(!hist) return 1.;\r
925 \r
926  //\r
927  TAxis *ax = hist->GetAxis(0);\r
928  TAxis *ay = hist->GetAxis(1);\r
929  TAxis *az = hist->GetAxis(2);\r
930 \r
931  Int_t binx = ax->FindBin(zv);\r
932  Int_t biny = ay->FindBin(pt);\r
933  Int_t binz = az->FindBin(eta);\r
934  Int_t dim[3] = {binx,biny,binz};\r
935 \r
936  Double_t fact  = hist->GetBinContent(dim);  \r
937 \r
938 return fact;\r
939 }\r
940 \r
941 //_____________________________________________________________________________\r
942 Double_t AlidNdPtCorrection::GetContFactZvPtEta(THnSparse *hist, Double_t zv, Double_t pt, Double_t eta) const {\r
943 // return contamination correction factor F(zv,pt,eta)\r
944 \r
945  if(!hist) return 1.0;\r
946 \r
947  //\r
948  TAxis *ax = hist->GetAxis(0);\r
949  TAxis *ay = hist->GetAxis(1);\r
950  TAxis *az = hist->GetAxis(2);\r
951 \r
952  Int_t binx = ax->FindBin(zv);\r
953  Int_t biny = ay->FindBin(pt);\r
954  Int_t binz = az->FindBin(eta);\r
955  Int_t dim[3] = {binx,biny,binz};\r
956 \r
957  Double_t fact  = 1.0-hist->GetBinContent(dim);  \r
958  //Double_t fact  = hist->GetBinContent(dim);  \r
959 \r
960 return fact;\r
961 }\r
962 \r
963 //_____________________________________________________________________________\r
964 Double_t AlidNdPtCorrection::GetCorrFactZvMult(THnSparse *hist, Double_t zv, Int_t mult) const {\r
965 // return correction factor F(zv,mult)\r
966 \r
967  if(!hist) return 1.;\r
968 \r
969  TAxis *ax = hist->GetAxis(0);\r
970  TAxis *ay = hist->GetAxis(1);\r
971  Int_t binx = ax->FindBin(zv);\r
972  Int_t biny = ay->FindBin(mult);\r
973  Int_t dim[2] = {binx,biny};\r
974 \r
975  Double_t fact  = hist->GetBinContent(dim);  \r
976 \r
977 \r
978 return fact;\r
979 }\r
980 \r
981 //_____________________________________________________________________________\r
982 Double_t AlidNdPtCorrection::GetContFactZvMult(THnSparse *hist, Double_t zv, Int_t mult) const {\r
983 // return contamination correction factor F(zv,mult)\r
984 \r
985  if(!hist) return 1.;\r
986 \r
987  TAxis *ax = hist->GetAxis(0);\r
988  TAxis *ay = hist->GetAxis(1);\r
989  Int_t binx = ax->FindBin(zv);\r
990  Int_t biny = ay->FindBin(mult);\r
991  Int_t dim[2] = {binx,biny};\r
992  Double_t fact  = 1.0-hist->GetBinContent(dim);  \r
993 \r
994 return fact;\r
995 }\r
996 \r
997 //_____________________________________________________________________________\r
998 Long64_t AlidNdPtCorrection::Merge(TCollection* list) \r
999 {\r
1000   // Merge list of objects (needed by PROOF)\r
1001 \r
1002   if (!list)\r
1003   return 0;\r
1004 \r
1005   if (list->IsEmpty())\r
1006   return 1;\r
1007 \r
1008   TIterator* iter = list->MakeIterator();\r
1009   TObject* obj = 0;\r
1010 \r
1011   // collection of generated histograms\r
1012 \r
1013   Int_t count=0;\r
1014   while((obj = iter->Next()) != 0) {\r
1015     AlidNdPtCorrection* entry = dynamic_cast<AlidNdPtCorrection*>(obj);\r
1016     if (entry == 0) continue; \r
1017   \r
1018     fMCEventHist1->Add(entry->fMCEventHist1);\r
1019     fRecEventHist1->Add(entry->fRecEventHist1);\r
1020 \r
1021     fMCAllEventMultHist1->Add(entry->fMCAllEventMultHist1);\r
1022     fMCAllNDEventMultHist1->Add(entry->fMCAllNDEventMultHist1);\r
1023     fMCAllNSDEventMultHist1->Add(entry->fMCAllNSDEventMultHist1);\r
1024     fMCTriggerMultHist1->Add(entry->fMCTriggerMultHist1);\r
1025     fMCEventMultHist1->Add(entry->fMCEventMultHist1);\r
1026 \r
1027     fMCAllPrimTrackMultHist1->Add(entry->fMCAllPrimTrackMultHist1);\r
1028     fMCNDEventAllPrimTrackMultHist1->Add(entry->fMCNDEventAllPrimTrackMultHist1);\r
1029     fMCNSDEventAllPrimTrackMultHist1->Add(entry->fMCNSDEventAllPrimTrackMultHist1);\r
1030     fMCTriggerPrimTrackMultHist1->Add(entry->fMCTriggerPrimTrackMultHist1);\r
1031     fMCEventPrimTrackMultHist1->Add(entry->fMCEventPrimTrackMultHist1);\r
1032 \r
1033     for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {\r
1034       fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);\r
1035     }\r
1036 \r
1037     for(Int_t i=0; i<8; i++) {\r
1038       fCorrRecTrackMultHist1[i]->Add(entry->fCorrRecTrackMultHist1[i]);\r
1039       fPtvsPt[i]->Add(entry->fPtvsPt[i]);\r
1040     }\r
1041 \r
1042     for(Int_t i=0; i<5; i++) {\r
1043       fCorrRecEventHist1[i]->Add(entry->fCorrRecEventHist1[i]);\r
1044       fCorrRecEventHist2[i]->Add(entry->fCorrRecEventHist2[i]);\r
1045     }\r
1046 \r
1047   count++;\r
1048   }\r
1049 \r
1050 return count;\r
1051 }\r
1052  \r
1053 Int_t AlidNdPtCorrection::GetTrueMult(THnSparse *hist, Int_t mult)\r
1054 {\r
1055  if(!hist) return 0;\r
1056  Int_t true_mult = 0;\r
1057 \r
1058  // 0 bins exluded\r
1059  TAxis *ax = hist->GetAxis(0);\r
1060  TAxis *ay = hist->GetAxis(1);\r
1061  ax->SetRange(1,ax->GetNbins());\r
1062  ay->SetRange(1,ay->GetNbins());\r
1063 \r
1064  // measured mult\r
1065  ax->SetRangeUser((Float_t)mult,(Float_t)mult); \r
1066 \r
1067  // get true multiplicity\r
1068  TH1D *h1 = (TH1D *)hist->Projection(1);\r
1069  true_mult = (Int_t)h1->GetMean();\r
1070 \r
1071  return true_mult;\r
1072 }\r
1073 \r
1074 //_____________________________________________________________________________\r
1075 void AlidNdPtCorrection::Analyse() \r
1076 {\r
1077   // Analyse histograms\r
1078   //\r
1079   TH1::AddDirectory(kFALSE);\r
1080   TH1 *h = 0, *hs=0, *hsc=0; \r
1081   TH2 *h2D = 0; \r
1082 \r
1083   TObjArray *aFolderObj = new TObjArray;\r
1084 \r
1085   //\r
1086   // get cuts\r
1087   //\r
1088   AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
1089   AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
1090   AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
1091 \r
1092   if(!evtCuts || !accCuts || !esdTrackCuts) {\r
1093     Error("AlidNdPtCutAnalysis::Analyse()", "cuts not available");\r
1094     return;\r
1095   }\r
1096 \r
1097   //\r
1098   // set min and max values\r
1099   //\r
1100   Double_t minPt = accCuts->GetMinPt();\r
1101   Double_t maxPt = accCuts->GetMaxPt();\r
1102   Double_t minEta = accCuts->GetMinEta();\r
1103   Double_t maxEta = accCuts->GetMaxEta()-0.00001;\r
1104 \r
1105   //\r
1106   // pt profile\r
1107   //\r
1108   char name[256];\r
1109   for(Int_t i=0; i<8; i++) {\r
1110     h2D = fPtvsPt[i]->Projection(1,0);\r
1111     sprintf(name,"PtvsMeanPt_%d",i);\r
1112     aFolderObj->Add(h2D);\r
1113   }\r
1114 \r
1115   //\r
1116   // event level \r
1117   //\r
1118   h = fCorrRecEventHist1[0]->Projection(1);\r
1119   h->SetName("mult_event_not_corrected");\r
1120   aFolderObj->Add(h);\r
1121 \r
1122   h = fCorrRecEventHist1[1]->Projection(1);\r
1123   h->SetName("mult_event_vertex_corrected");\r
1124   aFolderObj->Add(h);\r
1125 \r
1126   h = fCorrRecEventHist1[2]->Projection(1);\r
1127   h->SetName("mult_trigger_vertex_corrected");\r
1128   aFolderObj->Add(h);\r
1129 \r
1130   h = fCorrRecEventHist1[3]->Projection(1);\r
1131   h->SetName("mult_ND_trigger_vertex_corrected");\r
1132   aFolderObj->Add(h);\r
1133 \r
1134   h = fCorrRecEventHist1[4]->Projection(1);\r
1135   h->SetName("mult_NSD_trigger_vertex_corrected");\r
1136   aFolderObj->Add(h);\r
1137 \r
1138   // empty events\r
1139   h = fCorrRecEventHist2[0]->Projection(1);\r
1140   h->SetName("mult_empty_event_not_corrected");\r
1141   aFolderObj->Add(h);\r
1142 \r
1143   h = fCorrRecEventHist2[1]->Projection(1);\r
1144   h->SetName("mult_empty_event_vertex_corrected");\r
1145   aFolderObj->Add(h);\r
1146 \r
1147   h = fCorrRecEventHist2[2]->Projection(1);\r
1148   h->SetName("mult_empty_trigger_vertex_corrected");\r
1149   aFolderObj->Add(h);\r
1150 \r
1151   h = fCorrRecEventHist2[3]->Projection(1);\r
1152   h->SetName("mult_empty_ND_trigger_vertex_corrected");\r
1153   aFolderObj->Add(h);\r
1154 \r
1155   h = fCorrRecEventHist2[4]->Projection(1);\r
1156   h->SetName("mult_empty_NSD_trigger_vertex_corrected");\r
1157   aFolderObj->Add(h);\r
1158  \r
1159   //\r
1160   // MC available\r
1161   //\r
1162   if(IsUseMCInfo()) {\r
1163 \r
1164   // mc \r
1165   h = fMCAllEventMultHist1->Projection(1);\r
1166   h->SetName("mc_mult_event_acc_prim");\r
1167   aFolderObj->Add(h);\r
1168 \r
1169   h = fMCAllNDEventMultHist1->Projection(1);\r
1170   h->SetName("mc_mult_ND_event_acc_prim");\r
1171   aFolderObj->Add(h);\r
1172 \r
1173   h = fMCAllNSDEventMultHist1->Projection(1);\r
1174   h->SetName("mc_mult_NSD_event_acc_prim");\r
1175   aFolderObj->Add(h);\r
1176 \r
1177   h = fMCTriggerMultHist1->Projection(1);\r
1178   h->SetName("mc_mult_trigger_acc_prim");\r
1179   aFolderObj->Add(h);\r
1180 \r
1181   h = fMCEventMultHist1->Projection(1);\r
1182   h->SetName("mc_mult_trigger_event_acc_prim");\r
1183   aFolderObj->Add(h);\r
1184 \r
1185 \r
1186   //\r
1187   // track level\r
1188   //\r
1189   \r
1190   // limit eta range\r
1191   for(Int_t i=0;i<8;i++) { \r
1192       fCorrRecTrackMultHist1[i]->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1193       fCorrRecTrackMultHist1[i]->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
1194   }\r
1195   fMCAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1196   fMCAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
1197 \r
1198   fMCNDEventAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1199   fMCNDEventAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
1200 \r
1201   fMCNSDEventAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1202   fMCNSDEventAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
1203 \r
1204   fMCTriggerPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1205   fMCTriggerPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
1206 \r
1207   fMCEventPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1208   fMCEventPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
1209 \r
1210   } // end use MC info \r
1211   \r
1212   //\r
1213   h2D = fCorrRecTrackMultHist1[3]->Projection(1,0);\r
1214   h2D->SetName("pt_eta_rec_track_mult_eff_cont_corrected");\r
1215   aFolderObj->Add(h2D);\r
1216 \r
1217   h2D = fCorrRecTrackMultHist1[4]->Projection(1,0);\r
1218   h2D->SetName("pt_eta_rec_event_track_mult_eff_cont_corrected");\r
1219   aFolderObj->Add(h2D);\r
1220 \r
1221   h2D = fCorrRecTrackMultHist1[5]->Projection(1,0);\r
1222   h2D->SetName("pt_eta_rec_trig_event_track_mult_eff_cont_corrected");\r
1223   aFolderObj->Add(h2D);\r
1224 \r
1225   h2D = fCorrRecTrackMultHist1[6]->Projection(1,0);\r
1226   h2D->SetName("pt_eta_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1227   aFolderObj->Add(h2D);\r
1228 \r
1229   h2D = fCorrRecTrackMultHist1[7]->Projection(1,0);\r
1230   h2D->SetName("pt_eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1231   aFolderObj->Add(h2D);\r
1232 \r
1233 \r
1234   //\r
1235   h2D = fCorrRecTrackMultHist1[3]->Projection(2,0);\r
1236   h2D->SetName("pt_mult_rec_track_mult_eff_cont_corrected");\r
1237   aFolderObj->Add(h2D);\r
1238 \r
1239   h2D = fCorrRecTrackMultHist1[4]->Projection(2,0);\r
1240   h2D->SetName("pt_mult_rec_event_track_mult_eff_cont_corrected");\r
1241   aFolderObj->Add(h2D);\r
1242 \r
1243   h2D = fCorrRecTrackMultHist1[5]->Projection(2,0);\r
1244   h2D->SetName("pt_mult_rec_trig_event_track_mult_eff_cont_corrected");\r
1245   aFolderObj->Add(h2D);\r
1246 \r
1247   h2D = fCorrRecTrackMultHist1[6]->Projection(2,0);\r
1248   h2D->SetName("pt_mult_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1249   aFolderObj->Add(h2D);\r
1250 \r
1251   h2D = fCorrRecTrackMultHist1[7]->Projection(2,0);\r
1252   h2D->SetName("pt_mult_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1253   aFolderObj->Add(h2D);\r
1254 \r
1255   // pt axis\r
1256 \r
1257   h = fCorrRecTrackMultHist1[0]->Projection(0);\r
1258   h->SetName("pt_rec_track_not_corrected");\r
1259   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1260   hs->SetName("pt_rec_track_not_corrected_s");\r
1261   aFolderObj->Add(hs);\r
1262 \r
1263   //\r
1264   h = fCorrRecTrackMultHist1[1]->Projection(0);\r
1265   h->SetName("pt_rec_track_cont_corrected");\r
1266   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1267   hs->SetName("pt_rec_track_cont_corrected_s");\r
1268   aFolderObj->Add(hs);\r
1269 \r
1270   hsc = (TH1D*)hs->Clone("pt_rec_track_cont_corr_fact");\r
1271   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_not_corrected_s"));\r
1272   aFolderObj->Add(hsc);\r
1273 \r
1274   //\r
1275   h = fCorrRecTrackMultHist1[2]->Projection(0);\r
1276   h->SetName("pt_rec_track_eff_cont_corrected");\r
1277   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1278   hs->SetName("pt_rec_track_eff_cont_corrected_s");\r
1279   aFolderObj->Add(hs);\r
1280 \r
1281   hsc = (TH1D*)hs->Clone("pt_rec_track_eff_corr_fact");\r
1282   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_cont_corrected_s"));\r
1283   aFolderObj->Add(hsc);\r
1284 \r
1285   //\r
1286   h = fCorrRecTrackMultHist1[3]->Projection(0);\r
1287   h->SetName("pt_rec_track_mult_eff_cont_corrected");\r
1288   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1289   hs->SetName("pt_rec_track_mult_eff_cont_corrected_s");\r
1290   aFolderObj->Add(hs);\r
1291 \r
1292   hsc = (TH1D*)hs->Clone("pt_rec_track_mult_corr_fact");\r
1293   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_eff_cont_corrected_s"));\r
1294   aFolderObj->Add(hsc);\r
1295 \r
1296   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1297   hsc = (TH1D*)hs->Clone("pt_rec_track_mult_eff_cont_corr_fact");\r
1298   hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_not_corrected_s"));\r
1299   aFolderObj->Add(hsc);\r
1300 \r
1301   hsc = (TH1D*)hs->Clone();\r
1302   hsc->SetName("pt_rec_track_mult_eff_cont_corrected_s_norm");\r
1303   hsc->Scale(1./(fCorrRecEventHist1[0]->Projection(1)->Integral()));\r
1304   aFolderObj->Add(hsc);\r
1305 \r
1306   //\r
1307   h = fCorrRecTrackMultHist1[4]->Projection(0);\r
1308   h->SetName("pt_rec_event_track_mult_eff_cont_corrected");\r
1309   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1310   hs->SetName("pt_rec_event_track_mult_eff_cont_corrected_s");\r
1311   aFolderObj->Add(hs);\r
1312 \r
1313   hsc = (TH1D*)hs->Clone();\r
1314   hsc->SetName("pt_rec_event_track_mult_eff_cont_corrected_s_norm");\r
1315   hsc->Scale(1./(fCorrRecEventHist1[1]->Projection(1)->Integral()+fCorrRecEventHist2[1]->Projection(1)->Integral()));\r
1316   aFolderObj->Add(hsc);\r
1317 \r
1318   //\r
1319   h = fCorrRecTrackMultHist1[5]->Projection(0);\r
1320   h->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected");\r
1321   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1322   hs->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s");\r
1323   aFolderObj->Add(hs);\r
1324 \r
1325   hsc = (TH1D*)hs->Clone();\r
1326   hsc->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1327   hsc->Scale(1./(fCorrRecEventHist1[2]->Projection(1)->Integral() + fCorrRecEventHist2[2]->Projection(1)->Integral()));\r
1328   aFolderObj->Add(hsc);\r
1329 \r
1330   //\r
1331   h = fCorrRecTrackMultHist1[6]->Projection(0);\r
1332   h->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1333   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1334   hs->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s");\r
1335   aFolderObj->Add(hs);\r
1336 \r
1337   hsc = (TH1D*)hs->Clone();\r
1338   hsc->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1339   hsc->Scale(1./(fCorrRecEventHist1[3]->Projection(1)->Integral()+fCorrRecEventHist2[3]->Projection(1)->Integral()));\r
1340   aFolderObj->Add(hsc);\r
1341 \r
1342   //\r
1343   h = fCorrRecTrackMultHist1[7]->Projection(0);\r
1344   h->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1345   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1346   hs->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s");\r
1347   aFolderObj->Add(hs);\r
1348 \r
1349   hsc = (TH1D*)hs->Clone();\r
1350   hsc->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1351   hsc->Scale(1./(fCorrRecEventHist1[4]->Projection(1)->Integral()+fCorrRecEventHist2[4]->Projection(1)->Integral()));\r
1352   aFolderObj->Add(hsc);\r
1353 \r
1354   // eta axis\r
1355   h = fCorrRecTrackMultHist1[0]->Projection(1);\r
1356   h->SetName("eta_rec_track_not_corrected");\r
1357   aFolderObj->Add(h);\r
1358   \r
1359   h = fCorrRecTrackMultHist1[1]->Projection(1);\r
1360   h->SetName("eta_rec_track_cont_corrected");\r
1361   aFolderObj->Add(h);\r
1362 \r
1363   h = fCorrRecTrackMultHist1[2]->Projection(1);\r
1364   h->SetName("eta_rec_track_eff_cont_corrected");\r
1365   aFolderObj->Add(h);\r
1366 \r
1367   h = fCorrRecTrackMultHist1[3]->Projection(1);\r
1368   h->SetName("eta_rec_track_mult_eff_cont_corrected");\r
1369   aFolderObj->Add(h);\r
1370 \r
1371   h = fCorrRecTrackMultHist1[4]->Projection(1);\r
1372   h->SetName("eta_rec_event_track_mult_eff_cont_corrected");\r
1373   aFolderObj->Add(h);\r
1374 \r
1375   h = fCorrRecTrackMultHist1[5]->Projection(1);\r
1376   h->SetName("eta_rec_trig_event_track_mult_eff_cont_corrected");\r
1377   aFolderObj->Add(h);\r
1378 \r
1379   h = fCorrRecTrackMultHist1[6]->Projection(1);\r
1380   h->SetName("eta_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1381   aFolderObj->Add(h);\r
1382 \r
1383   h = fCorrRecTrackMultHist1[7]->Projection(1);\r
1384   h->SetName("eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1385   aFolderObj->Add(h);\r
1386 \r
1387 \r
1388   //\r
1389   // MC available\r
1390   //\r
1391   if(IsUseMCInfo()) {\r
1392 \r
1393   //\r
1394   h2D = fMCAllPrimTrackMultHist1->Projection(2,0);\r
1395   h2D->SetName("mc_all_pt_mult_acc_prim");\r
1396   aFolderObj->Add(h2D);\r
1397 \r
1398   h2D = fMCAllPrimTrackMultHist1->Projection(1,0);\r
1399   h2D->SetName("mc_all_eta_pt_acc_prim");\r
1400   aFolderObj->Add(h2D);\r
1401 \r
1402   h2D = fMCNDEventAllPrimTrackMultHist1->Projection(2,0);\r
1403   h2D->SetName("mc_ND_all_pt_mult_acc_prim");\r
1404   aFolderObj->Add(h2D);\r
1405 \r
1406   h2D = fMCNDEventAllPrimTrackMultHist1->Projection(1,0);\r
1407   h2D->SetName("mc_ND_all_eta_pt_acc_prim");\r
1408   aFolderObj->Add(h2D);\r
1409 \r
1410   h2D = fMCNSDEventAllPrimTrackMultHist1->Projection(2,0);\r
1411   h2D->SetName("mc_NSD_all_pt_mult_acc_prim");\r
1412   aFolderObj->Add(h2D);\r
1413 \r
1414   h2D = fMCNSDEventAllPrimTrackMultHist1->Projection(1,0);\r
1415   h2D->SetName("mc_NSD_all_eta_pt_acc_prim");\r
1416   aFolderObj->Add(h2D);\r
1417 \r
1418   //\r
1419 \r
1420   h = fMCAllPrimTrackMultHist1->Projection(0);\r
1421   h->SetName("mc_all_pt_acc_prim");\r
1422   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1423   hs->SetName("mc_all_pt_acc_prim_s");\r
1424   aFolderObj->Add(hs);\r
1425 \r
1426   hsc = (TH1D*)hs->Clone();\r
1427   hsc->SetName("mc_all_pt_acc_prim_s_norm");\r
1428   hsc->Scale(1./fMCAllEventMultHist1->Projection(1)->Integral());\r
1429   aFolderObj->Add(hsc);\r
1430 \r
1431   h = fMCNDEventAllPrimTrackMultHist1->Projection(0);\r
1432   h->SetName("mc_ND_all_pt_acc_prim");\r
1433   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1434   hs->SetName("mc_ND_all_pt_acc_prim_s");\r
1435   aFolderObj->Add(hs);\r
1436 \r
1437   hsc = (TH1D*)hs->Clone();\r
1438   hsc->SetName("mc_ND_all_pt_acc_prim_s_norm");\r
1439   hsc->Scale(1./fMCAllNDEventMultHist1->Projection(1)->Integral());\r
1440   aFolderObj->Add(hsc);\r
1441 \r
1442   h = fMCNSDEventAllPrimTrackMultHist1->Projection(0);\r
1443   h->SetName("mc_NSD_all_pt_acc_prim");\r
1444   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1445   hs->SetName("mc_NSD_all_pt_acc_prim_s");\r
1446   aFolderObj->Add(hs);\r
1447 \r
1448   hsc = (TH1D*)hs->Clone();\r
1449   hsc->SetName("mc_NSD_all_pt_acc_prim_s_norm");\r
1450   hsc->Scale(1./fMCAllNSDEventMultHist1->Projection(1)->Integral());\r
1451   aFolderObj->Add(hsc);\r
1452 \r
1453   h = fMCTriggerPrimTrackMultHist1->Projection(0);\r
1454   h->SetName("mc_trigger_all_pt_acc_prim");\r
1455   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1456   hs->SetName("mc_trigger_all_pt_acc_prim_s");\r
1457   aFolderObj->Add(hs);\r
1458 \r
1459   hsc = (TH1D*)hs->Clone();\r
1460   hsc->SetName("mc_trigger_all_pt_acc_prim_s_norm");\r
1461   hsc->Scale(1./fMCTriggerMultHist1->Projection(1)->Integral());\r
1462   aFolderObj->Add(hsc);\r
1463 \r
1464   h = fMCEventPrimTrackMultHist1->Projection(0);\r
1465   h->SetName("mc_all_pt_acc_trig_event_prim");\r
1466   hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1467   hs->SetName("mc_all_pt_acc_trig_event_prim_s");\r
1468   aFolderObj->Add(hs);\r
1469 \r
1470   hsc = (TH1D*)hs->Clone();\r
1471   hsc->SetName("mc_all_pt_acc_trig_event_prim_s_norm");\r
1472   hsc->Scale(1./fMCEventMultHist1->Projection(1)->Integral());\r
1473   aFolderObj->Add(hsc);\r
1474 \r
1475   //\r
1476 \r
1477   h = fMCAllPrimTrackMultHist1->Projection(1);\r
1478   h->SetName("mc_all_eta_acc_prim");\r
1479   aFolderObj->Add(h);\r
1480 \r
1481   h = fMCNDEventAllPrimTrackMultHist1->Projection(1);\r
1482   h->SetName("mc_ND_all_eta_acc_prim");\r
1483   aFolderObj->Add(h);\r
1484 \r
1485   h = fMCNSDEventAllPrimTrackMultHist1->Projection(1);\r
1486   h->SetName("mc_NSD_all_eta_acc_prim");\r
1487   aFolderObj->Add(h);\r
1488 \r
1489   h = fMCTriggerPrimTrackMultHist1->Projection(1);\r
1490   h->SetName("mc_trigger_all_eta_acc_prim");\r
1491   aFolderObj->Add(h);\r
1492 \r
1493   h = fMCEventPrimTrackMultHist1->Projection(1);\r
1494   h->SetName("mc_all_eta_acc_trig_event_prim");\r
1495   aFolderObj->Add(h);\r
1496 \r
1497   //\r
1498   // calculate ratios (rec / mc)\r
1499   //\r
1500   \r
1501   hs = (TH1*)aFolderObj->FindObject("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1502   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_trig_event_track_mult_eff_cont_corrected"); \r
1503   hsc->Sumw2();\r
1504   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_pt_acc_prim_s_norm"));\r
1505   aFolderObj->Add(hsc);\r
1506 \r
1507   hs = (TH1*)aFolderObj->FindObject("eta_rec_trig_event_track_mult_eff_cont_corrected");\r
1508   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_trig_event_track_mult_eff_cont_corrected"); \r
1509   hsc->Sumw2();\r
1510   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_eta_acc_prim"));\r
1511   aFolderObj->Add(hsc);\r
1512 \r
1513   //\r
1514   hs = (TH1*)aFolderObj->FindObject("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1515   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_ND_trig_event_track_mult_eff_cont_corrected"); \r
1516   hsc->Sumw2();\r
1517   hsc->Divide((TH1*)aFolderObj->FindObject("mc_ND_all_pt_acc_prim_s_norm"));\r
1518   aFolderObj->Add(hsc);\r
1519 \r
1520   hs = (TH1*)aFolderObj->FindObject("eta_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1521   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_ND_trig_event_track_mult_eff_cont_corrected"); \r
1522   hsc->Sumw2();\r
1523   hsc->Divide((TH1*)aFolderObj->FindObject("mc_ND_all_eta_acc_prim"));\r
1524   aFolderObj->Add(hsc);\r
1525 \r
1526   //\r
1527   hs = (TH1*)aFolderObj->FindObject("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1528   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_NSD_trig_event_track_mult_eff_cont_corrected"); \r
1529   hsc->Sumw2();\r
1530   hsc->Divide((TH1*)aFolderObj->FindObject("mc_NSD_all_pt_acc_prim_s_norm"));\r
1531   aFolderObj->Add(hsc);\r
1532 \r
1533   hs = (TH1*)aFolderObj->FindObject("eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1534   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_NSD_trig_event_track_mult_eff_cont_corrected"); \r
1535   hsc->Sumw2();\r
1536   hsc->Divide((TH1*)aFolderObj->FindObject("mc_NSD_all_eta_acc_prim"));\r
1537   aFolderObj->Add(hsc);\r
1538 \r
1539   //\r
1540   hs = (TH1*)aFolderObj->FindObject("pt_rec_event_track_mult_eff_cont_corrected_s_norm");\r
1541   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_event_track_mult_eff_cont_corrected"); \r
1542   hsc->Sumw2();\r
1543   hsc->Divide((TH1*)aFolderObj->FindObject("mc_trigger_all_pt_acc_prim_s_norm"));\r
1544   aFolderObj->Add(hsc);\r
1545 \r
1546   hs = (TH1*)aFolderObj->FindObject("eta_rec_event_track_mult_eff_cont_corrected");\r
1547   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_event_track_mult_eff_cont_corrected"); \r
1548   hsc->Sumw2();\r
1549   hsc->Divide((TH1*)aFolderObj->FindObject("mc_trigger_all_eta_acc_prim"));\r
1550   aFolderObj->Add(hsc);\r
1551 \r
1552   // track level\r
1553   hs = (TH1*)aFolderObj->FindObject("pt_rec_track_mult_eff_cont_corrected_s_norm");\r
1554   hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_track_mult_eff_cont_corrected"); \r
1555   hsc->Sumw2();\r
1556   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_pt_acc_trig_event_prim_s_norm"));\r
1557   aFolderObj->Add(hsc);\r
1558 \r
1559   hs = (TH1*)aFolderObj->FindObject("eta_rec_track_mult_eff_cont_corrected");\r
1560   hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_track_mult_eff_cont_corrected"); \r
1561   hsc->Sumw2();\r
1562   hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_eta_acc_trig_event_prim"));\r
1563   aFolderObj->Add(hsc);\r
1564 \r
1565   } // end MC infor available\r
1566 \r
1567   // export objects to analysis folder\r
1568   fCorrectionFolder = ExportToFolder(aFolderObj);\r
1569 \r
1570   // delete only TObjArray\r
1571   if(aFolderObj) delete aFolderObj;\r
1572 }\r
1573 \r
1574 //_____________________________________________________________________________\r
1575 TFolder* AlidNdPtCorrection::ExportToFolder(TObjArray * array) \r
1576 {\r
1577   // recreate folder avery time and export objects to new one\r
1578   //\r
1579   AlidNdPtCorrection * comp=this;\r
1580   TFolder *folder = comp->GetCorrectionFolder();\r
1581 \r
1582   TString name, title;\r
1583   TFolder *newFolder = 0;\r
1584   Int_t i = 0;\r
1585   Int_t size = array->GetSize();\r
1586 \r
1587   if(folder) { \r
1588      // get name and title from old folder\r
1589      name = folder->GetName();  \r
1590      title = folder->GetTitle();  \r
1591 \r
1592          // delete old one\r
1593      delete folder;\r
1594 \r
1595          // create new one\r
1596      newFolder = CreateFolder(name.Data(),title.Data());\r
1597      newFolder->SetOwner();\r
1598 \r
1599          // add objects to folder\r
1600      while(i < size) {\r
1601            newFolder->Add(array->At(i));\r
1602            i++;\r
1603          }\r
1604   }\r
1605 \r
1606 return newFolder;\r
1607 }\r
1608 \r
1609 //_____________________________________________________________________________\r
1610 TFolder* AlidNdPtCorrection::CreateFolder(TString name,TString title) { \r
1611 // create folder for analysed histograms\r
1612 //\r
1613 TFolder *folder = 0;\r
1614   folder = new TFolder(name.Data(),title.Data());\r
1615 \r
1616   return folder;\r
1617 }\r