classes for limnosity related info from LHC DIP data.
[u/mrichter/AliRoot.git] / PWG0 / dNdPt / AlidNdPtCorrection.cxx
CommitLineData
0aaa8b91 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
bad4ba69 34#include "AliPWG0Helper.h"\r
35#include "AlidNdPtHelper.h"\r
0aaa8b91 36#include "AlidNdPtCorrection.h"\r
37\r
38using namespace std;\r
39\r
40ClassImp(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
0aaa8b91 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
92AlidNdPtCorrection::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
0aaa8b91 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
bad4ba69 129 fPtvsPt[i] = 0;\r
0aaa8b91 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
141AlidNdPtCorrection::~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
bad4ba69 166 if(fPtvsPt[i]) delete fPtvsPt[i]; fPtvsPt[i]=0;\r
0aaa8b91 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
0aaa8b91 175}\r
176\r
177//_____________________________________________________________________________\r
178void 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
0aaa8b91 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
bad4ba69 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
0aaa8b91 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
bad4ba69 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
0aaa8b91 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
505void 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
bad4ba69 515 // get selection cuts\r
0aaa8b91 516 AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
517 AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
518 AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
519\r
bad4ba69 520 if(!evtCuts || !accCuts || !esdTrackCuts) {\r
0aaa8b91 521 AliDebug(AliLog::kError, "cuts not available");\r
522 return;\r
523 }\r
524\r
bad4ba69 525 // trigger selection\r
526 Bool_t isEventTriggered = kTRUE;\r
527 if(evtCuts->IsTriggerRequired()) {\r
f537848a 528 if(IsUseMCInfo()) { \r
529 static AliTriggerAnalysis* triggerAnalysis = new AliTriggerAnalysis;\r
530 isEventTriggered = triggerAnalysis->IsTriggerFired(esdEvent, GetTrigger());\r
531 }\r
532 else {\r
533 isEventTriggered = esdEvent->IsTriggerClassFired(GetTriggerClass());\r
534 }\r
bad4ba69 535 }\r
0aaa8b91 536\r
bad4ba69 537 // use MC information\r
0aaa8b91 538 AliHeader* header = 0;\r
539 AliGenEventHeader* genHeader = 0;\r
540 AliStack* stack = 0;\r
541 TArrayF vtxMC(3);\r
bad4ba69 542 AliPWG0Helper::MCProcessType evtType = AliPWG0Helper::kInvalidProcess;\r
0aaa8b91 543 Int_t multMCTrueTracks = 0;\r
544\r
545 if(IsUseMCInfo())\r
546 {\r
547 if(!mcEvent) {\r
548 AliDebug(AliLog::kError, "mcEvent not available");\r
549 return;\r
550 }\r
551\r
552 // get MC event header\r
553 header = mcEvent->Header();\r
554 if (!header) {\r
555 AliDebug(AliLog::kError, "Header not available");\r
556 return;\r
557 }\r
558 // MC particle stack\r
559 stack = mcEvent->Stack();\r
560 if (!stack) {\r
561 AliDebug(AliLog::kError, "Stack not available");\r
562 return;\r
563 }\r
564\r
565 // get event type (ND=0x1, DD=0x2, SD=0x4)\r
bad4ba69 566 evtType = AliPWG0Helper::GetEventProcessType(header);\r
0aaa8b91 567 //Printf("evtType %d \n", evtType);\r
568 AliDebug(AliLog::kDebug+1, Form("Found process type %d", evtType));\r
569\r
570 // get MC vertex\r
571 genHeader = header->GenEventHeader();\r
572 if (!genHeader) {\r
573 AliDebug(AliLog::kError, "Could not retrieve genHeader from Header");\r
574 return;\r
575 }\r
576 genHeader->PrimaryVertex(vtxMC);\r
577\r
578 // Fill MC event histogram\r
579 Double_t vMCEventHist1[3]={vtxMC[0],vtxMC[1],vtxMC[2]};\r
580 fMCEventHist1->Fill(vMCEventHist1);\r
581\r
582 // multipliticy of all MC primary tracks\r
583 // in Zvtx, pt and eta ranges\r
584 multMCTrueTracks = AlidNdPtHelper::GetMCTrueTrackMult(mcEvent,evtCuts,accCuts);\r
585\r
586 } // end bUseMC\r
587\r
588 // get reconstructed vertex \r
589 const AliESDVertex* vtxESD = 0; \r
590 Bool_t isRecVertex = kFALSE;\r
591 if(evtCuts->IsRecVertexRequired()) \r
592 {\r
bad4ba69 593 Bool_t bRedoTPCVertex = evtCuts->IsRedoTPCVertex();\r
594 Bool_t bUseConstraints = evtCuts->IsUseBeamSpotConstraint();\r
595 vtxESD = AlidNdPtHelper::GetVertex(esdEvent,evtCuts,accCuts,esdTrackCuts,GetAnalysisMode(),kFALSE,bRedoTPCVertex,bUseConstraints); \r
847e74b2 596 isRecVertex = AlidNdPtHelper::TestRecVertex(vtxESD, GetAnalysisMode(), kFALSE); \r
0aaa8b91 597 }\r
598 if( IsUseMCInfo() && !evtCuts->IsRecVertexRequired() ) {\r
599 vtxESD = new AliESDVertex(vtxMC[2],10.,genHeader->NProduced());\r
600 isRecVertex = kTRUE;\r
601 }\r
602 Bool_t isEventOK = evtCuts->AcceptEvent(esdEvent,mcEvent,vtxESD) && isRecVertex; \r
603 //printf("isEventOK %d \n",isEventOK);\r
604\r
605 //\r
847e74b2 606 // get multiplicity vertex contributors\r
0aaa8b91 607 //\r
608 Int_t multMBTracks = 0; \r
609 if(GetAnalysisMode() == AlidNdPtHelper::kTPC || GetAnalysisMode() == AlidNdPtHelper::kMCRec) { \r
610 multMBTracks = AlidNdPtHelper::GetTPCMBTrackMult(esdEvent,evtCuts,accCuts,esdTrackCuts);\r
611 } \r
847e74b2 612 else if(GetAnalysisMode() == AlidNdPtHelper::kTPCSPDvtx || GetAnalysisMode() == AlidNdPtHelper::kTPCSPDvtxUpdate || GetAnalysisMode() == AlidNdPtHelper::kMCRec) {\r
613 //multMBTracks = AlidNdPtHelper::GetSPDMBTrackMult(esdEvent,0.0);\r
614 multMBTracks = vtxESD->GetNContributors();\r
0aaa8b91 615 } \r
616 else {\r
617 AliDebug(AliLog::kError, Form("Found analysis type %s", GetAnalysisMode()));\r
618 return; \r
619 }\r
620\r
621 //\r
622 // correct event and track histograms\r
623 //\r
624 TObjArray *allChargedTracks=0;\r
625 Int_t multRec=0, multRecTemp=0;\r
626 Int_t *labelsRec=0;\r
627\r
628 if(isEventOK && isEventTriggered)\r
629 {\r
630 // get all charged tracks\r
bad4ba69 631 //allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,vtxESD,GetAnalysisMode());\r
632 allChargedTracks = AlidNdPtHelper::GetAllChargedTracks(esdEvent,GetAnalysisMode());\r
0aaa8b91 633 if(!allChargedTracks) return;\r
634\r
635 Int_t entries = allChargedTracks->GetEntries();\r
636 labelsRec = new Int_t[entries];\r
637\r
638 // calculate mult of reconstructed tracks\r
639 for(Int_t i=0; i<entries;++i) \r
640 {\r
641 AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
642 if(!track) continue;\r
847e74b2 643 if(track->Charge()==0) continue;\r
644\r
645 // only postive charged \r
646 if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0) \r
647 continue;\r
648 \r
649 // only negative charged \r
650 if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0) \r
651 continue;\r
652\r
653 if(esdTrackCuts->AcceptTrack(track) && accCuts->AcceptTrack(track)) \r
654 multRecTemp++;\r
0aaa8b91 655 } \r
656\r
847e74b2 657 //\r
0aaa8b91 658 for(Int_t i=0; i<entries;++i) \r
659 {\r
660 AliESDtrack *track = (AliESDtrack*)allChargedTracks->At(i);\r
661 if(!track) continue;\r
847e74b2 662 if(track->Charge()==0) continue;\r
663\r
664 // only postive charged \r
665 if(GetParticleMode() == AlidNdPtHelper::kPlus && track->Charge() < 0) \r
666 continue;\r
667 \r
668 // only negative charged \r
669 if(GetParticleMode() == AlidNdPtHelper::kMinus && track->Charge() > 0) \r
670 continue;\r
0aaa8b91 671 \r
847e74b2 672 // track-level corrections\r
673 if(esdTrackCuts->AcceptTrack(track)) \r
0aaa8b91 674 {\r
847e74b2 675 if (GetAnalysisMode()==AlidNdPtHelper::kTPCSPDvtxUpdate) \r
676 {\r
677 //\r
678 // update track parameters\r
679 //\r
680 AliExternalTrackParam cParam;\r
681 track->RelateToVertexTPC(esdEvent->GetPrimaryVertexSPD(),esdEvent->GetMagneticField(),kVeryBig,&cParam);\r
682 track->Set(cParam.GetX(),cParam.GetAlpha(),cParam.GetParameter(),cParam.GetCovariance());\r
683\r
684 if(accCuts->AcceptTrack(track)) { \r
685 FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxESD->GetZv(),multRecTemp); \r
686 labelsRec[multRec] = TMath::Abs(track->GetLabel());\r
687 multRec++;\r
688 }\r
689 }\r
690 else if (GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) \r
691 { \r
692 //\r
693 // Replace rec with MC\r
694 //\r
695 if(accCuts->AcceptTrack(track)) { \r
696 FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxMC[2],multRecTemp); \r
697 labelsRec[multRec] = TMath::Abs(track->GetLabel());\r
698 multRec++;\r
699 }\r
700 } \r
701 else {\r
702 //\r
703 // all the rest tracking scenarios\r
704 //\r
705 if(accCuts->AcceptTrack(track)) { \r
706 FillHistograms(track,stack,AlidNdPtHelper::kRecTracks,vtxESD->GetZv(),multRecTemp); \r
707 labelsRec[multRec] = TMath::Abs(track->GetLabel());\r
708 multRec++;\r
709 }\r
710 }\r
711 }\r
0aaa8b91 712 }\r
713 // event-level corrections\r
714 if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) { \r
715 FillHistograms(AlidNdPtHelper::kRecEvents,vtxMC[2],multMBTracks);\r
716 }\r
717 else {\r
718 FillHistograms(AlidNdPtHelper::kRecEvents,vtxESD->GetZv(),multMBTracks);\r
719 }\r
720\r
721 // control event histograms\r
722 Double_t vRecEventHist1[3] = {vtxESD->GetXv(),vtxESD->GetYv(),vtxESD->GetZv()};\r
723 fRecEventHist1->Fill(vRecEventHist1);\r
724 } \r
725\r
726 // empty events corrections\r
727 // no reconstructed zv\r
728 if(isEventTriggered && multMBTracks==0) {\r
729 if(GetAnalysisMode()==AlidNdPtHelper::kMCRec && IsUseMCInfo()) {\r
730 FillHistograms(AlidNdPtHelper::kTriggeredEvents,vtxMC[2],multMBTracks);\r
731 }\r
732 else {\r
733 Double_t zv = fZvNorm->GetRandom();\r
734 FillHistograms(AlidNdPtHelper::kTriggeredEvents,zv,multMBTracks);\r
735 }\r
736 }\r
737\r
738 if(IsUseMCInfo()) \r
739 {\r
740 // select MC events \r
741 if(evtCuts->AcceptMCEvent(mcEvent))\r
742 {\r
743 //\r
744 // event histograms\r
745 //\r
746 Double_t vMCEventMatrix[2] = {vtxMC[2],multMBTracks};\r
747 fMCAllEventMultHist1->Fill(vMCEventMatrix);\r
bad4ba69 748 if(evtType == AliPWG0Helper::kND) {\r
0aaa8b91 749 fMCAllNDEventMultHist1->Fill(vMCEventMatrix);\r
750 }\r
bad4ba69 751 if(evtType != AliPWG0Helper::kSD) {\r
0aaa8b91 752 fMCAllNSDEventMultHist1->Fill(vMCEventMatrix);\r
753 }\r
754 if(isEventTriggered) fMCTriggerMultHist1->Fill(vMCEventMatrix);\r
755 if(isEventTriggered && isEventOK) fMCEventMultHist1->Fill(vMCEventMatrix);\r
756\r
757 //\r
758 // MC histograms for efficiency studies\r
759 //\r
760 Int_t nPart = stack->GetNtrack();\r
761 for (Int_t iMc = 0; iMc < nPart; ++iMc) \r
762 {\r
763 // print MC stack info\r
764 //AlidNdPtHelper::PrintMCInfo(stack,iMc);\r
765\r
766 TParticle* particle = stack->Particle(iMc);\r
767 if (!particle)\r
768 continue;\r
769\r
770 // only charged particles\r
771 Double_t charge = particle->GetPDG()->Charge()/3.;\r
772 if (charge == 0.0)\r
773 continue;\r
847e74b2 774\r
775 // only postive charged \r
776 if(GetParticleMode() == AlidNdPtHelper::kPlus && charge < 0.) \r
777 continue;\r
778 \r
779 // only negative charged \r
780 if(GetParticleMode() == AlidNdPtHelper::kMinus && charge > 0.) \r
781 continue;\r
0aaa8b91 782 \r
783 // physical primary\r
784 Bool_t prim = stack->IsPhysicalPrimary(iMc);\r
785\r
786 // all primaries in acceptance\r
787 if(!accCuts->AcceptTrack(particle)) continue;\r
788 if(!prim) continue;\r
789\r
790 Double_t gpt = particle->Pt();\r
791 Double_t geta = particle->Eta();\r
792\r
793 Double_t valMCAllTrackMultHist1[3] = {gpt,geta,multRec}; \r
794 fMCAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
bad4ba69 795 if(evtType == AliPWG0Helper::kND) {\r
0aaa8b91 796 fMCNDEventAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
797 }\r
bad4ba69 798 if(evtType != AliPWG0Helper::kSD) {\r
0aaa8b91 799 fMCNSDEventAllPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
800 }\r
801 if(isEventTriggered) fMCTriggerPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
802 if(isEventTriggered && isEventOK) fMCEventPrimTrackMultHist1->Fill(valMCAllTrackMultHist1);\r
803 }\r
804 }\r
805 } // end bUseMC\r
806\r
807 if(allChargedTracks) delete allChargedTracks; allChargedTracks = 0;\r
808 if(labelsRec) delete [] labelsRec; labelsRec = 0;\r
809\r
f537848a 810 if(!evtCuts->IsRecVertexRequired() && vtxESD != NULL) delete vtxESD;\r
0aaa8b91 811}\r
812\r
813//_____________________________________________________________________________\r
814void AlidNdPtCorrection::FillHistograms(AlidNdPtHelper::EventObject eventObj, Double_t zv, Int_t multMBTracks)\r
815{\r
816 //\r
817 // Fill corrected histograms\r
818 //\r
819\r
820 Double_t vEventMatrix[2] = {zv,multMBTracks};\r
821 //\r
822 // Correct for efficiency \r
823 //\r
824 if(eventObj == AlidNdPtHelper::kRecEvents && multMBTracks>0) \r
825 {\r
826 Double_t corrToMBF = GetCorrFactZvMult(fCorrEventMatrix,zv,multMBTracks);\r
827 Double_t corrToInelF = GetCorrFactZvMult(fCorrTriggerMBtoInelEventMatrix,zv,multMBTracks);\r
828 Double_t corrToNDF = GetCorrFactZvMult(fCorrTriggerMBtoNDEventMatrix,zv,multMBTracks);\r
829 Double_t corrToNSDF = GetCorrFactZvMult(fCorrTriggerMBtoNSDEventMatrix,zv,multMBTracks);\r
830 //printf("corrToMBF %f, corrToInelF %f, corrToNDF %f corrToNSDF %f \n",corrToMBF,corrToInelF,corrToNDF,corrToNSDF);\r
831\r
832 fCorrRecEventHist1[0]->Fill(vEventMatrix);\r
833 fCorrRecEventHist1[1]->Fill(vEventMatrix,corrToMBF);\r
834 fCorrRecEventHist1[2]->Fill(vEventMatrix,corrToMBF*corrToInelF);\r
835 fCorrRecEventHist1[3]->Fill(vEventMatrix,corrToMBF*corrToNDF);\r
836 fCorrRecEventHist1[4]->Fill(vEventMatrix,corrToMBF*corrToNSDF);\r
837 }\r
838\r
839 if(eventObj==AlidNdPtHelper::kTriggeredEvents && multMBTracks==0) // empty triggered events\r
840 {\r
841 Int_t bin = fZvEmptyEventsNorm->FindBin(zv);\r
842 Double_t Fz = fZvEmptyEventsNorm->GetBinContent(bin);\r
843 Double_t corrToInelF0 = GetCorrFactZvMult(fCorrTriggerMBtoInelEventMatrix,zv,multMBTracks);\r
844 Double_t corrToNDF0 = GetCorrFactZvMult(fCorrTriggerMBtoNDEventMatrix,zv,multMBTracks);\r
845 Double_t corrToNSDF0 = GetCorrFactZvMult(fCorrTriggerMBtoNSDEventMatrix,zv,multMBTracks);\r
846 //printf("Fz %f, corrToInelF0 %f, corrToNDF0 %f, corrToNSDF0 %f \n",Fz,corrToInelF0,corrToNDF0,corrToNSDF0);\r
847\r
848 fCorrRecEventHist2[0]->Fill(vEventMatrix);\r
849 fCorrRecEventHist2[1]->Fill(vEventMatrix,Fz);\r
850 fCorrRecEventHist2[2]->Fill(vEventMatrix,Fz*corrToInelF0);\r
851 fCorrRecEventHist2[3]->Fill(vEventMatrix,Fz*corrToNDF0);\r
852 fCorrRecEventHist2[4]->Fill(vEventMatrix,Fz*corrToNSDF0);\r
853 }\r
854}\r
855\r
856//_____________________________________________________________________________\r
857void AlidNdPtCorrection::FillHistograms(AliESDtrack *esdTrack, AliStack *stack, AlidNdPtHelper::TrackObject trackObj, Double_t zv, Int_t mult)\r
858{\r
859 //\r
860 // Fill ESD track and MC histograms \r
861 //\r
862 if(!esdTrack) return;\r
863\r
864 //Float_t q = esdTrack->Charge();\r
865 Float_t pt = esdTrack->Pt();\r
866 Float_t eta = esdTrack->Eta();\r
867 Float_t phi = esdTrack->Phi();\r
868\r
869 if(stack && GetAnalysisMode() == AlidNdPtHelper::kMCRec) \r
870 {\r
871 Int_t label = TMath::Abs(esdTrack->GetLabel());\r
872 \r
873 TParticle* particle = stack->Particle(label);\r
874 if(!particle) return;\r
875 \r
876 Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3\r
877 if(gq==0) return;\r
878 Float_t gpt = particle->Pt();\r
879 Float_t geta = particle->Eta();\r
880 Float_t gphi = particle->Phi();\r
881\r
882 // replace reconstructed values with MC\r
883 pt = gpt;\r
884 eta = geta;\r
885 phi = gphi;\r
886 }\r
887\r
888 //\r
889 // Fill histograms\r
890 //\r
891 Double_t values[3] = {pt,eta,phi}; \r
892 fRecTrackHist1[trackObj]->Fill(values);\r
893\r
894 //\r
895 // Correct for contamination and efficiency \r
896 //\r
897 if(trackObj == AlidNdPtHelper::kRecTracks || GetAnalysisMode() == AlidNdPtHelper::kMCRec) \r
898 {\r
899 // track level corrections\r
900 Double_t trackEffF = GetCorrFactZvPtEta(fCorrTrackMatrix,zv,pt,eta);\r
901 Double_t trackContF = GetContFactZvPtEta(fContTrackMatrix,zv,pt,eta);\r
902 Double_t multTrackContF = GetContFactZvPtEta(fContMultTrackMatrix,zv,pt,eta);\r
903 //printf("zv %f, pt %f, eta %f \n",zv,pt,eta);\r
904 //printf("trackEffF %f, trackContF %f, multTrackContF %f \n", trackEffF, trackContF, multTrackContF);\r
905 \r
906 // track-event level corrections\r
907 Double_t vertexEffF = GetCorrFactZvPtEta(fCorrTrackEventMatrix,zv,pt,eta);\r
908 Double_t trigMBToInel = GetCorrFactZvPtEta(fCorrTriggerMBtoInelTrackEventMatrix,zv,pt,eta); \r
909 Double_t trigMBToND = GetCorrFactZvPtEta(fCorrTriggerMBtoNDTrackEventMatrix,zv,pt,eta);\r
910 Double_t trigMBToNSD = GetCorrFactZvPtEta(fCorrTriggerMBtoNSDTrackEventMatrix,zv,pt,eta);\r
911 //printf("vertexEffF %f, trigMBToInel %f, trigMBToNSD %f \n", vertexEffF, trigMBToInel, trigMBToNSD);\r
bad4ba69 912 \r
913 Double_t corrF[8] = { 1.0, \r
914 trackContF,\r
915 trackContF*trackEffF,\r
916 trackContF*trackEffF*multTrackContF,\r
917 trackContF*trackEffF*multTrackContF*vertexEffF,\r
918 trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToInel,\r
919 trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToND,\r
920 trackContF*trackEffF*multTrackContF*vertexEffF*trigMBToNSD\r
921 }; \r
922 \r
923 // Fill histograms\r
0aaa8b91 924 Double_t valCorrRecTrackMultHist1[3] = {pt,eta,mult}; \r
bad4ba69 925 Double_t valPtvsPt[2] = {pt,pt}; \r
926 for(Int_t i=0; i<8; i++) {\r
927 fCorrRecTrackMultHist1[i]->Fill(valCorrRecTrackMultHist1,corrF[i]);\r
928 fPtvsPt[i]->Fill(valPtvsPt,corrF[i]);\r
929 }\r
0aaa8b91 930 }\r
931}\r
932\r
933//_____________________________________________________________________________\r
934void AlidNdPtCorrection::FillHistograms(AliStack *stack, Int_t /*label*/, AlidNdPtHelper::TrackObject /*trackObj*/, Int_t /*mult*/)\r
935{\r
936 // Fill MC histograms\r
937 if(!stack) return;\r
938\r
939 /*\r
940 TParticle* particle = stack->Particle(label);\r
941 if(!particle) return;\r
942\r
943 Int_t mother_pdg = -1;\r
944 TParticle* mother = 0;\r
945\r
946 //TParticle* prim_mother = AlidNdPtHelper::FindPrimaryMother(stack,label);\r
947 Int_t motherLabel = particle->GetMother(0); \r
948 if(motherLabel>0) mother = stack->Particle(motherLabel);\r
949 if(mother) mother_pdg = TMath::Abs(mother->GetPdgCode()); // take abs for visualisation only\r
950 Int_t mech = particle->GetUniqueID(); // production mechanism\r
951\r
952 Double_t gq = particle->GetPDG()->Charge()/3.0; // Charge units |e|/3 \r
953 Float_t gpt = particle->Pt();\r
954 Float_t qgpt = particle->Pt() * gq;\r
955 Float_t geta = particle->Eta();\r
956 Float_t gphi = particle->Phi();\r
957 Float_t gpz = particle->Pz();\r
958\r
959 Bool_t prim = stack->IsPhysicalPrimary(label);\r
960 Float_t vx = particle->Vx(); Float_t vy = particle->Vy(); Float_t vz = particle->Vz();\r
961\r
962 Int_t pid=-1;\r
963 if (TMath::Abs(particle->GetPdgCode()) == kElectron) { pid = 0; }\r
964 else if (TMath::Abs(particle->GetPdgCode()) == kMuonMinus) { pid = 1; }\r
965 else if (TMath::Abs(particle->GetPdgCode()) == kPiPlus) { pid = 2; }\r
966 else if (TMath::Abs(particle->GetPdgCode()) == kKPlus) { pid = 3; }\r
967 else if (TMath::Abs(particle->GetPdgCode()) == kProton) { pid = 4; }\r
968 else { pid = 5; }\r
969 */\r
970\r
971 //if(!prim) printf("prim_mother %d, mother %d, particle %d, production mech %d\n",prim_mother->GetPdgCode(),mother->GetPdgCode(), particle->GetPdgCode(),mech);\r
972 \r
973}\r
974\r
975//_____________________________________________________________________________\r
976Double_t AlidNdPtCorrection::GetCorrFactZvPtEta(THnSparse *hist, Double_t zv, Double_t pt, Double_t eta) const {\r
977// return correction factor F(zv,pt,eta)\r
978\r
979 if(!hist) return 1.;\r
980\r
981 //\r
982 TAxis *ax = hist->GetAxis(0);\r
983 TAxis *ay = hist->GetAxis(1);\r
984 TAxis *az = hist->GetAxis(2);\r
985\r
986 Int_t binx = ax->FindBin(zv);\r
987 Int_t biny = ay->FindBin(pt);\r
988 Int_t binz = az->FindBin(eta);\r
989 Int_t dim[3] = {binx,biny,binz};\r
990\r
991 Double_t fact = hist->GetBinContent(dim); \r
992\r
993return fact;\r
994}\r
995\r
996//_____________________________________________________________________________\r
997Double_t AlidNdPtCorrection::GetContFactZvPtEta(THnSparse *hist, Double_t zv, Double_t pt, Double_t eta) const {\r
998// return contamination correction factor F(zv,pt,eta)\r
999\r
1000 if(!hist) return 1.0;\r
1001\r
1002 //\r
1003 TAxis *ax = hist->GetAxis(0);\r
1004 TAxis *ay = hist->GetAxis(1);\r
1005 TAxis *az = hist->GetAxis(2);\r
1006\r
1007 Int_t binx = ax->FindBin(zv);\r
1008 Int_t biny = ay->FindBin(pt);\r
1009 Int_t binz = az->FindBin(eta);\r
1010 Int_t dim[3] = {binx,biny,binz};\r
1011\r
1012 Double_t fact = 1.0-hist->GetBinContent(dim); \r
1013 //Double_t fact = hist->GetBinContent(dim); \r
1014\r
1015return fact;\r
1016}\r
1017\r
1018//_____________________________________________________________________________\r
1019Double_t AlidNdPtCorrection::GetCorrFactZvMult(THnSparse *hist, Double_t zv, Int_t mult) const {\r
1020// return correction factor F(zv,mult)\r
1021\r
1022 if(!hist) return 1.;\r
1023\r
1024 TAxis *ax = hist->GetAxis(0);\r
1025 TAxis *ay = hist->GetAxis(1);\r
1026 Int_t binx = ax->FindBin(zv);\r
1027 Int_t biny = ay->FindBin(mult);\r
1028 Int_t dim[2] = {binx,biny};\r
1029\r
1030 Double_t fact = hist->GetBinContent(dim); \r
1031\r
1032\r
1033return fact;\r
1034}\r
1035\r
1036//_____________________________________________________________________________\r
1037Double_t AlidNdPtCorrection::GetContFactZvMult(THnSparse *hist, Double_t zv, Int_t mult) const {\r
1038// return contamination correction factor F(zv,mult)\r
1039\r
1040 if(!hist) return 1.;\r
1041\r
1042 TAxis *ax = hist->GetAxis(0);\r
1043 TAxis *ay = hist->GetAxis(1);\r
1044 Int_t binx = ax->FindBin(zv);\r
1045 Int_t biny = ay->FindBin(mult);\r
1046 Int_t dim[2] = {binx,biny};\r
1047 Double_t fact = 1.0-hist->GetBinContent(dim); \r
1048\r
1049return fact;\r
1050}\r
1051\r
1052//_____________________________________________________________________________\r
1053Long64_t AlidNdPtCorrection::Merge(TCollection* list) \r
1054{\r
1055 // Merge list of objects (needed by PROOF)\r
1056\r
1057 if (!list)\r
1058 return 0;\r
1059\r
1060 if (list->IsEmpty())\r
1061 return 1;\r
1062\r
1063 TIterator* iter = list->MakeIterator();\r
1064 TObject* obj = 0;\r
1065\r
1066 // collection of generated histograms\r
1067\r
1068 Int_t count=0;\r
1069 while((obj = iter->Next()) != 0) {\r
1070 AlidNdPtCorrection* entry = dynamic_cast<AlidNdPtCorrection*>(obj);\r
1071 if (entry == 0) continue; \r
1072 \r
1073 fMCEventHist1->Add(entry->fMCEventHist1);\r
1074 fRecEventHist1->Add(entry->fRecEventHist1);\r
1075\r
1076 fMCAllEventMultHist1->Add(entry->fMCAllEventMultHist1);\r
1077 fMCAllNDEventMultHist1->Add(entry->fMCAllNDEventMultHist1);\r
1078 fMCAllNSDEventMultHist1->Add(entry->fMCAllNSDEventMultHist1);\r
1079 fMCTriggerMultHist1->Add(entry->fMCTriggerMultHist1);\r
1080 fMCEventMultHist1->Add(entry->fMCEventMultHist1);\r
1081\r
1082 fMCAllPrimTrackMultHist1->Add(entry->fMCAllPrimTrackMultHist1);\r
1083 fMCNDEventAllPrimTrackMultHist1->Add(entry->fMCNDEventAllPrimTrackMultHist1);\r
1084 fMCNSDEventAllPrimTrackMultHist1->Add(entry->fMCNSDEventAllPrimTrackMultHist1);\r
1085 fMCTriggerPrimTrackMultHist1->Add(entry->fMCTriggerPrimTrackMultHist1);\r
1086 fMCEventPrimTrackMultHist1->Add(entry->fMCEventPrimTrackMultHist1);\r
1087\r
1088 for(Int_t i=0; i<AlidNdPtHelper::kCutSteps; i++) {\r
1089 fRecTrackHist1[i]->Add(entry->fRecTrackHist1[i]);\r
1090 }\r
1091\r
1092 for(Int_t i=0; i<8; i++) {\r
1093 fCorrRecTrackMultHist1[i]->Add(entry->fCorrRecTrackMultHist1[i]);\r
bad4ba69 1094 fPtvsPt[i]->Add(entry->fPtvsPt[i]);\r
0aaa8b91 1095 }\r
1096\r
1097 for(Int_t i=0; i<5; i++) {\r
1098 fCorrRecEventHist1[i]->Add(entry->fCorrRecEventHist1[i]);\r
1099 fCorrRecEventHist2[i]->Add(entry->fCorrRecEventHist2[i]);\r
1100 }\r
1101\r
0aaa8b91 1102 count++;\r
1103 }\r
1104\r
1105return count;\r
1106}\r
1107 \r
1108Int_t AlidNdPtCorrection::GetTrueMult(THnSparse *hist, Int_t mult)\r
1109{\r
1110 if(!hist) return 0;\r
1111 Int_t true_mult = 0;\r
1112\r
1113 // 0 bins exluded\r
1114 TAxis *ax = hist->GetAxis(0);\r
1115 TAxis *ay = hist->GetAxis(1);\r
1116 ax->SetRange(1,ax->GetNbins());\r
1117 ay->SetRange(1,ay->GetNbins());\r
1118\r
1119 // measured mult\r
1120 ax->SetRangeUser((Float_t)mult,(Float_t)mult); \r
1121\r
1122 // get true multiplicity\r
1123 TH1D *h1 = (TH1D *)hist->Projection(1);\r
1124 true_mult = (Int_t)h1->GetMean();\r
1125\r
1126 return true_mult;\r
1127}\r
1128\r
1129//_____________________________________________________________________________\r
1130void AlidNdPtCorrection::Analyse() \r
1131{\r
1132 // Analyse histograms\r
1133 //\r
1134 TH1::AddDirectory(kFALSE);\r
1135 TH1 *h = 0, *hs=0, *hsc=0; \r
1136 TH2 *h2D = 0; \r
1137\r
1138 TObjArray *aFolderObj = new TObjArray;\r
1139\r
1140 //\r
1141 // get cuts\r
1142 //\r
1143 AlidNdPtEventCuts *evtCuts = GetEventCuts(); \r
1144 AlidNdPtAcceptanceCuts *accCuts = GetAcceptanceCuts(); \r
1145 AliESDtrackCuts *esdTrackCuts = GetTrackCuts(); \r
1146\r
1147 if(!evtCuts || !accCuts || !esdTrackCuts) {\r
1148 Error("AlidNdPtCutAnalysis::Analyse()", "cuts not available");\r
1149 return;\r
1150 }\r
1151\r
1152 //\r
1153 // set min and max values\r
1154 //\r
847e74b2 1155 //Double_t minPt = accCuts->GetMinPt();\r
1156 //Double_t maxPt = accCuts->GetMaxPt();\r
1157 //Double_t minEta = accCuts->GetMinEta();\r
1158 //Double_t maxEta = accCuts->GetMaxEta()-0.00001;\r
0aaa8b91 1159\r
1160 //\r
1161 // pt profile\r
1162 //\r
bad4ba69 1163 char name[256];\r
1164 for(Int_t i=0; i<8; i++) {\r
1165 h2D = fPtvsPt[i]->Projection(1,0);\r
1166 sprintf(name,"PtvsMeanPt_%d",i);\r
1167 aFolderObj->Add(h2D);\r
1168 }\r
0aaa8b91 1169\r
1170 //\r
1171 // event level \r
1172 //\r
1173 h = fCorrRecEventHist1[0]->Projection(1);\r
1174 h->SetName("mult_event_not_corrected");\r
1175 aFolderObj->Add(h);\r
1176\r
1177 h = fCorrRecEventHist1[1]->Projection(1);\r
1178 h->SetName("mult_event_vertex_corrected");\r
1179 aFolderObj->Add(h);\r
1180\r
1181 h = fCorrRecEventHist1[2]->Projection(1);\r
1182 h->SetName("mult_trigger_vertex_corrected");\r
1183 aFolderObj->Add(h);\r
1184\r
1185 h = fCorrRecEventHist1[3]->Projection(1);\r
1186 h->SetName("mult_ND_trigger_vertex_corrected");\r
1187 aFolderObj->Add(h);\r
1188\r
1189 h = fCorrRecEventHist1[4]->Projection(1);\r
1190 h->SetName("mult_NSD_trigger_vertex_corrected");\r
1191 aFolderObj->Add(h);\r
1192\r
1193 // empty events\r
1194 h = fCorrRecEventHist2[0]->Projection(1);\r
1195 h->SetName("mult_empty_event_not_corrected");\r
1196 aFolderObj->Add(h);\r
1197\r
1198 h = fCorrRecEventHist2[1]->Projection(1);\r
1199 h->SetName("mult_empty_event_vertex_corrected");\r
1200 aFolderObj->Add(h);\r
1201\r
1202 h = fCorrRecEventHist2[2]->Projection(1);\r
1203 h->SetName("mult_empty_trigger_vertex_corrected");\r
1204 aFolderObj->Add(h);\r
1205\r
1206 h = fCorrRecEventHist2[3]->Projection(1);\r
1207 h->SetName("mult_empty_ND_trigger_vertex_corrected");\r
1208 aFolderObj->Add(h);\r
1209\r
1210 h = fCorrRecEventHist2[4]->Projection(1);\r
1211 h->SetName("mult_empty_NSD_trigger_vertex_corrected");\r
1212 aFolderObj->Add(h);\r
bad4ba69 1213 \r
1214 //\r
1215 // MC available\r
1216 //\r
1217 if(IsUseMCInfo()) {\r
0aaa8b91 1218\r
1219 // mc \r
1220 h = fMCAllEventMultHist1->Projection(1);\r
1221 h->SetName("mc_mult_event_acc_prim");\r
1222 aFolderObj->Add(h);\r
1223\r
1224 h = fMCAllNDEventMultHist1->Projection(1);\r
1225 h->SetName("mc_mult_ND_event_acc_prim");\r
1226 aFolderObj->Add(h);\r
1227\r
1228 h = fMCAllNSDEventMultHist1->Projection(1);\r
1229 h->SetName("mc_mult_NSD_event_acc_prim");\r
1230 aFolderObj->Add(h);\r
1231\r
1232 h = fMCTriggerMultHist1->Projection(1);\r
1233 h->SetName("mc_mult_trigger_acc_prim");\r
1234 aFolderObj->Add(h);\r
1235\r
1236 h = fMCEventMultHist1->Projection(1);\r
1237 h->SetName("mc_mult_trigger_event_acc_prim");\r
1238 aFolderObj->Add(h);\r
1239\r
1240\r
1241 //\r
1242 // track level\r
1243 //\r
1244 \r
1245 // limit eta range\r
1246 for(Int_t i=0;i<8;i++) { \r
3405c7c2 1247 //fCorrRecTrackMultHist1[i]->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1248 //fCorrRecTrackMultHist1[i]->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
0aaa8b91 1249 }\r
3405c7c2 1250 //fMCAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1251 //fMCAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
0aaa8b91 1252\r
3405c7c2 1253 //fMCNDEventAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1254 //fMCNDEventAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
0aaa8b91 1255\r
3405c7c2 1256 //fMCNSDEventAllPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1257 //fMCNSDEventAllPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
0aaa8b91 1258\r
3405c7c2 1259 //fMCTriggerPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1260 //fMCTriggerPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
0aaa8b91 1261\r
3405c7c2 1262 //fMCEventPrimTrackMultHist1->GetAxis(0)->SetRangeUser(minPt,maxPt);\r
1263 //fMCEventPrimTrackMultHist1->GetAxis(1)->SetRangeUser(minEta,maxEta);\r
bad4ba69 1264\r
1265 } // end use MC info \r
0aaa8b91 1266 \r
1267 //\r
1268 h2D = fCorrRecTrackMultHist1[3]->Projection(1,0);\r
1269 h2D->SetName("pt_eta_rec_track_mult_eff_cont_corrected");\r
1270 aFolderObj->Add(h2D);\r
1271\r
1272 h2D = fCorrRecTrackMultHist1[4]->Projection(1,0);\r
1273 h2D->SetName("pt_eta_rec_event_track_mult_eff_cont_corrected");\r
1274 aFolderObj->Add(h2D);\r
1275\r
1276 h2D = fCorrRecTrackMultHist1[5]->Projection(1,0);\r
1277 h2D->SetName("pt_eta_rec_trig_event_track_mult_eff_cont_corrected");\r
1278 aFolderObj->Add(h2D);\r
1279\r
1280 h2D = fCorrRecTrackMultHist1[6]->Projection(1,0);\r
1281 h2D->SetName("pt_eta_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1282 aFolderObj->Add(h2D);\r
1283\r
1284 h2D = fCorrRecTrackMultHist1[7]->Projection(1,0);\r
1285 h2D->SetName("pt_eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1286 aFolderObj->Add(h2D);\r
1287\r
1288\r
1289 //\r
1290 h2D = fCorrRecTrackMultHist1[3]->Projection(2,0);\r
1291 h2D->SetName("pt_mult_rec_track_mult_eff_cont_corrected");\r
1292 aFolderObj->Add(h2D);\r
1293\r
1294 h2D = fCorrRecTrackMultHist1[4]->Projection(2,0);\r
1295 h2D->SetName("pt_mult_rec_event_track_mult_eff_cont_corrected");\r
1296 aFolderObj->Add(h2D);\r
1297\r
1298 h2D = fCorrRecTrackMultHist1[5]->Projection(2,0);\r
1299 h2D->SetName("pt_mult_rec_trig_event_track_mult_eff_cont_corrected");\r
1300 aFolderObj->Add(h2D);\r
1301\r
1302 h2D = fCorrRecTrackMultHist1[6]->Projection(2,0);\r
1303 h2D->SetName("pt_mult_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1304 aFolderObj->Add(h2D);\r
1305\r
1306 h2D = fCorrRecTrackMultHist1[7]->Projection(2,0);\r
1307 h2D->SetName("pt_mult_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1308 aFolderObj->Add(h2D);\r
1309\r
1310 // pt axis\r
1311\r
1312 h = fCorrRecTrackMultHist1[0]->Projection(0);\r
1313 h->SetName("pt_rec_track_not_corrected");\r
1314 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1315 hs->SetName("pt_rec_track_not_corrected_s");\r
1316 aFolderObj->Add(hs);\r
1317\r
1318 //\r
1319 h = fCorrRecTrackMultHist1[1]->Projection(0);\r
1320 h->SetName("pt_rec_track_cont_corrected");\r
1321 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1322 hs->SetName("pt_rec_track_cont_corrected_s");\r
1323 aFolderObj->Add(hs);\r
1324\r
1325 hsc = (TH1D*)hs->Clone("pt_rec_track_cont_corr_fact");\r
1326 hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_not_corrected_s"));\r
1327 aFolderObj->Add(hsc);\r
1328\r
1329 //\r
1330 h = fCorrRecTrackMultHist1[2]->Projection(0);\r
1331 h->SetName("pt_rec_track_eff_cont_corrected");\r
1332 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1333 hs->SetName("pt_rec_track_eff_cont_corrected_s");\r
1334 aFolderObj->Add(hs);\r
1335\r
1336 hsc = (TH1D*)hs->Clone("pt_rec_track_eff_corr_fact");\r
1337 hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_cont_corrected_s"));\r
1338 aFolderObj->Add(hsc);\r
1339\r
1340 //\r
1341 h = fCorrRecTrackMultHist1[3]->Projection(0);\r
1342 h->SetName("pt_rec_track_mult_eff_cont_corrected");\r
1343 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1344 hs->SetName("pt_rec_track_mult_eff_cont_corrected_s");\r
1345 aFolderObj->Add(hs);\r
1346\r
1347 hsc = (TH1D*)hs->Clone("pt_rec_track_mult_corr_fact");\r
1348 hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_eff_cont_corrected_s"));\r
1349 aFolderObj->Add(hsc);\r
1350\r
1351 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1352 hsc = (TH1D*)hs->Clone("pt_rec_track_mult_eff_cont_corr_fact");\r
1353 hsc->Divide((TH1D*)aFolderObj->FindObject("pt_rec_track_not_corrected_s"));\r
1354 aFolderObj->Add(hsc);\r
1355\r
1356 hsc = (TH1D*)hs->Clone();\r
1357 hsc->SetName("pt_rec_track_mult_eff_cont_corrected_s_norm");\r
1358 hsc->Scale(1./(fCorrRecEventHist1[0]->Projection(1)->Integral()));\r
1359 aFolderObj->Add(hsc);\r
1360\r
1361 //\r
1362 h = fCorrRecTrackMultHist1[4]->Projection(0);\r
1363 h->SetName("pt_rec_event_track_mult_eff_cont_corrected");\r
1364 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1365 hs->SetName("pt_rec_event_track_mult_eff_cont_corrected_s");\r
1366 aFolderObj->Add(hs);\r
1367\r
1368 hsc = (TH1D*)hs->Clone();\r
1369 hsc->SetName("pt_rec_event_track_mult_eff_cont_corrected_s_norm");\r
1370 hsc->Scale(1./(fCorrRecEventHist1[1]->Projection(1)->Integral()+fCorrRecEventHist2[1]->Projection(1)->Integral()));\r
1371 aFolderObj->Add(hsc);\r
1372\r
1373 //\r
1374 h = fCorrRecTrackMultHist1[5]->Projection(0);\r
1375 h->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected");\r
1376 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1377 hs->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s");\r
1378 aFolderObj->Add(hs);\r
1379\r
1380 hsc = (TH1D*)hs->Clone();\r
1381 hsc->SetName("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1382 hsc->Scale(1./(fCorrRecEventHist1[2]->Projection(1)->Integral() + fCorrRecEventHist2[2]->Projection(1)->Integral()));\r
1383 aFolderObj->Add(hsc);\r
1384\r
1385 //\r
1386 h = fCorrRecTrackMultHist1[6]->Projection(0);\r
1387 h->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1388 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1389 hs->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s");\r
1390 aFolderObj->Add(hs);\r
1391\r
1392 hsc = (TH1D*)hs->Clone();\r
1393 hsc->SetName("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1394 hsc->Scale(1./(fCorrRecEventHist1[3]->Projection(1)->Integral()+fCorrRecEventHist2[3]->Projection(1)->Integral()));\r
1395 aFolderObj->Add(hsc);\r
1396\r
1397 //\r
1398 h = fCorrRecTrackMultHist1[7]->Projection(0);\r
1399 h->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1400 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1401 hs->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s");\r
1402 aFolderObj->Add(hs);\r
1403\r
1404 hsc = (TH1D*)hs->Clone();\r
1405 hsc->SetName("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1406 hsc->Scale(1./(fCorrRecEventHist1[4]->Projection(1)->Integral()+fCorrRecEventHist2[4]->Projection(1)->Integral()));\r
1407 aFolderObj->Add(hsc);\r
1408\r
1409 // eta axis\r
1410 h = fCorrRecTrackMultHist1[0]->Projection(1);\r
1411 h->SetName("eta_rec_track_not_corrected");\r
1412 aFolderObj->Add(h);\r
1413 \r
1414 h = fCorrRecTrackMultHist1[1]->Projection(1);\r
1415 h->SetName("eta_rec_track_cont_corrected");\r
1416 aFolderObj->Add(h);\r
1417\r
1418 h = fCorrRecTrackMultHist1[2]->Projection(1);\r
1419 h->SetName("eta_rec_track_eff_cont_corrected");\r
1420 aFolderObj->Add(h);\r
1421\r
1422 h = fCorrRecTrackMultHist1[3]->Projection(1);\r
1423 h->SetName("eta_rec_track_mult_eff_cont_corrected");\r
1424 aFolderObj->Add(h);\r
1425\r
1426 h = fCorrRecTrackMultHist1[4]->Projection(1);\r
1427 h->SetName("eta_rec_event_track_mult_eff_cont_corrected");\r
1428 aFolderObj->Add(h);\r
1429\r
1430 h = fCorrRecTrackMultHist1[5]->Projection(1);\r
1431 h->SetName("eta_rec_trig_event_track_mult_eff_cont_corrected");\r
1432 aFolderObj->Add(h);\r
1433\r
1434 h = fCorrRecTrackMultHist1[6]->Projection(1);\r
1435 h->SetName("eta_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1436 aFolderObj->Add(h);\r
1437\r
1438 h = fCorrRecTrackMultHist1[7]->Projection(1);\r
1439 h->SetName("eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1440 aFolderObj->Add(h);\r
1441\r
1442\r
bad4ba69 1443 //\r
1444 // MC available\r
1445 //\r
1446 if(IsUseMCInfo()) {\r
0aaa8b91 1447\r
1448 //\r
1449 h2D = fMCAllPrimTrackMultHist1->Projection(2,0);\r
1450 h2D->SetName("mc_all_pt_mult_acc_prim");\r
1451 aFolderObj->Add(h2D);\r
1452\r
1453 h2D = fMCAllPrimTrackMultHist1->Projection(1,0);\r
1454 h2D->SetName("mc_all_eta_pt_acc_prim");\r
1455 aFolderObj->Add(h2D);\r
1456\r
1457 h2D = fMCNDEventAllPrimTrackMultHist1->Projection(2,0);\r
1458 h2D->SetName("mc_ND_all_pt_mult_acc_prim");\r
1459 aFolderObj->Add(h2D);\r
1460\r
1461 h2D = fMCNDEventAllPrimTrackMultHist1->Projection(1,0);\r
1462 h2D->SetName("mc_ND_all_eta_pt_acc_prim");\r
1463 aFolderObj->Add(h2D);\r
1464\r
1465 h2D = fMCNSDEventAllPrimTrackMultHist1->Projection(2,0);\r
1466 h2D->SetName("mc_NSD_all_pt_mult_acc_prim");\r
1467 aFolderObj->Add(h2D);\r
1468\r
1469 h2D = fMCNSDEventAllPrimTrackMultHist1->Projection(1,0);\r
1470 h2D->SetName("mc_NSD_all_eta_pt_acc_prim");\r
1471 aFolderObj->Add(h2D);\r
1472\r
1473 //\r
1474\r
1475 h = fMCAllPrimTrackMultHist1->Projection(0);\r
1476 h->SetName("mc_all_pt_acc_prim");\r
1477 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1478 hs->SetName("mc_all_pt_acc_prim_s");\r
1479 aFolderObj->Add(hs);\r
1480\r
1481 hsc = (TH1D*)hs->Clone();\r
1482 hsc->SetName("mc_all_pt_acc_prim_s_norm");\r
1483 hsc->Scale(1./fMCAllEventMultHist1->Projection(1)->Integral());\r
1484 aFolderObj->Add(hsc);\r
1485\r
1486 h = fMCNDEventAllPrimTrackMultHist1->Projection(0);\r
1487 h->SetName("mc_ND_all_pt_acc_prim");\r
1488 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1489 hs->SetName("mc_ND_all_pt_acc_prim_s");\r
1490 aFolderObj->Add(hs);\r
1491\r
1492 hsc = (TH1D*)hs->Clone();\r
1493 hsc->SetName("mc_ND_all_pt_acc_prim_s_norm");\r
1494 hsc->Scale(1./fMCAllNDEventMultHist1->Projection(1)->Integral());\r
1495 aFolderObj->Add(hsc);\r
1496\r
1497 h = fMCNSDEventAllPrimTrackMultHist1->Projection(0);\r
1498 h->SetName("mc_NSD_all_pt_acc_prim");\r
1499 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1500 hs->SetName("mc_NSD_all_pt_acc_prim_s");\r
1501 aFolderObj->Add(hs);\r
1502\r
1503 hsc = (TH1D*)hs->Clone();\r
1504 hsc->SetName("mc_NSD_all_pt_acc_prim_s_norm");\r
1505 hsc->Scale(1./fMCAllNSDEventMultHist1->Projection(1)->Integral());\r
1506 aFolderObj->Add(hsc);\r
1507\r
1508 h = fMCTriggerPrimTrackMultHist1->Projection(0);\r
1509 h->SetName("mc_trigger_all_pt_acc_prim");\r
1510 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1511 hs->SetName("mc_trigger_all_pt_acc_prim_s");\r
1512 aFolderObj->Add(hs);\r
1513\r
1514 hsc = (TH1D*)hs->Clone();\r
1515 hsc->SetName("mc_trigger_all_pt_acc_prim_s_norm");\r
1516 hsc->Scale(1./fMCTriggerMultHist1->Projection(1)->Integral());\r
1517 aFolderObj->Add(hsc);\r
1518\r
1519 h = fMCEventPrimTrackMultHist1->Projection(0);\r
1520 h->SetName("mc_all_pt_acc_trig_event_prim");\r
1521 hs = AlidNdPtHelper::ScaleByBinWidth(h);\r
1522 hs->SetName("mc_all_pt_acc_trig_event_prim_s");\r
1523 aFolderObj->Add(hs);\r
1524\r
1525 hsc = (TH1D*)hs->Clone();\r
1526 hsc->SetName("mc_all_pt_acc_trig_event_prim_s_norm");\r
1527 hsc->Scale(1./fMCEventMultHist1->Projection(1)->Integral());\r
1528 aFolderObj->Add(hsc);\r
1529\r
1530 //\r
1531\r
1532 h = fMCAllPrimTrackMultHist1->Projection(1);\r
1533 h->SetName("mc_all_eta_acc_prim");\r
1534 aFolderObj->Add(h);\r
1535\r
1536 h = fMCNDEventAllPrimTrackMultHist1->Projection(1);\r
1537 h->SetName("mc_ND_all_eta_acc_prim");\r
1538 aFolderObj->Add(h);\r
1539\r
1540 h = fMCNSDEventAllPrimTrackMultHist1->Projection(1);\r
1541 h->SetName("mc_NSD_all_eta_acc_prim");\r
1542 aFolderObj->Add(h);\r
1543\r
1544 h = fMCTriggerPrimTrackMultHist1->Projection(1);\r
1545 h->SetName("mc_trigger_all_eta_acc_prim");\r
1546 aFolderObj->Add(h);\r
1547\r
1548 h = fMCEventPrimTrackMultHist1->Projection(1);\r
1549 h->SetName("mc_all_eta_acc_trig_event_prim");\r
1550 aFolderObj->Add(h);\r
1551\r
1552 //\r
1553 // calculate ratios (rec / mc)\r
1554 //\r
1555 \r
1556 hs = (TH1*)aFolderObj->FindObject("pt_rec_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1557 hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_trig_event_track_mult_eff_cont_corrected"); \r
1558 hsc->Sumw2();\r
1559 hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_pt_acc_prim_s_norm"));\r
1560 aFolderObj->Add(hsc);\r
1561\r
1562 hs = (TH1*)aFolderObj->FindObject("eta_rec_trig_event_track_mult_eff_cont_corrected");\r
1563 hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_trig_event_track_mult_eff_cont_corrected"); \r
1564 hsc->Sumw2();\r
1565 hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_eta_acc_prim"));\r
1566 aFolderObj->Add(hsc);\r
1567\r
1568 //\r
1569 hs = (TH1*)aFolderObj->FindObject("pt_rec_ND_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1570 hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_ND_trig_event_track_mult_eff_cont_corrected"); \r
1571 hsc->Sumw2();\r
1572 hsc->Divide((TH1*)aFolderObj->FindObject("mc_ND_all_pt_acc_prim_s_norm"));\r
1573 aFolderObj->Add(hsc);\r
1574\r
1575 hs = (TH1*)aFolderObj->FindObject("eta_rec_ND_trig_event_track_mult_eff_cont_corrected");\r
1576 hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_ND_trig_event_track_mult_eff_cont_corrected"); \r
1577 hsc->Sumw2();\r
1578 hsc->Divide((TH1*)aFolderObj->FindObject("mc_ND_all_eta_acc_prim"));\r
1579 aFolderObj->Add(hsc);\r
1580\r
1581 //\r
1582 hs = (TH1*)aFolderObj->FindObject("pt_rec_NSD_trig_event_track_mult_eff_cont_corrected_s_norm");\r
1583 hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_NSD_trig_event_track_mult_eff_cont_corrected"); \r
1584 hsc->Sumw2();\r
1585 hsc->Divide((TH1*)aFolderObj->FindObject("mc_NSD_all_pt_acc_prim_s_norm"));\r
1586 aFolderObj->Add(hsc);\r
1587\r
1588 hs = (TH1*)aFolderObj->FindObject("eta_rec_NSD_trig_event_track_mult_eff_cont_corrected");\r
1589 hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_NSD_trig_event_track_mult_eff_cont_corrected"); \r
1590 hsc->Sumw2();\r
1591 hsc->Divide((TH1*)aFolderObj->FindObject("mc_NSD_all_eta_acc_prim"));\r
1592 aFolderObj->Add(hsc);\r
1593\r
1594 //\r
1595 hs = (TH1*)aFolderObj->FindObject("pt_rec_event_track_mult_eff_cont_corrected_s_norm");\r
1596 hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_event_track_mult_eff_cont_corrected"); \r
1597 hsc->Sumw2();\r
1598 hsc->Divide((TH1*)aFolderObj->FindObject("mc_trigger_all_pt_acc_prim_s_norm"));\r
1599 aFolderObj->Add(hsc);\r
1600\r
1601 hs = (TH1*)aFolderObj->FindObject("eta_rec_event_track_mult_eff_cont_corrected");\r
1602 hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_event_track_mult_eff_cont_corrected"); \r
1603 hsc->Sumw2();\r
1604 hsc->Divide((TH1*)aFolderObj->FindObject("mc_trigger_all_eta_acc_prim"));\r
1605 aFolderObj->Add(hsc);\r
1606\r
1607 // track level\r
1608 hs = (TH1*)aFolderObj->FindObject("pt_rec_track_mult_eff_cont_corrected_s_norm");\r
1609 hsc = (TH1D*)hs->Clone("ratio_pt_rec_to_mc_track_mult_eff_cont_corrected"); \r
1610 hsc->Sumw2();\r
1611 hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_pt_acc_trig_event_prim_s_norm"));\r
1612 aFolderObj->Add(hsc);\r
1613\r
1614 hs = (TH1*)aFolderObj->FindObject("eta_rec_track_mult_eff_cont_corrected");\r
1615 hsc = (TH1D*)hs->Clone("ratio_eta_rec_to_mc_track_mult_eff_cont_corrected"); \r
1616 hsc->Sumw2();\r
1617 hsc->Divide((TH1*)aFolderObj->FindObject("mc_all_eta_acc_trig_event_prim"));\r
1618 aFolderObj->Add(hsc);\r
1619\r
bad4ba69 1620 } // end MC infor available\r
1621\r
0aaa8b91 1622 // export objects to analysis folder\r
1623 fCorrectionFolder = ExportToFolder(aFolderObj);\r
1624\r
1625 // delete only TObjArray\r
1626 if(aFolderObj) delete aFolderObj;\r
1627}\r
1628\r
1629//_____________________________________________________________________________\r
1630TFolder* AlidNdPtCorrection::ExportToFolder(TObjArray * array) \r
1631{\r
1632 // recreate folder avery time and export objects to new one\r
1633 //\r
1634 AlidNdPtCorrection * comp=this;\r
1635 TFolder *folder = comp->GetCorrectionFolder();\r
1636\r
1637 TString name, title;\r
1638 TFolder *newFolder = 0;\r
1639 Int_t i = 0;\r
1640 Int_t size = array->GetSize();\r
1641\r
1642 if(folder) { \r
1643 // get name and title from old folder\r
1644 name = folder->GetName(); \r
1645 title = folder->GetTitle(); \r
1646\r
1647 // delete old one\r
1648 delete folder;\r
1649\r
1650 // create new one\r
1651 newFolder = CreateFolder(name.Data(),title.Data());\r
1652 newFolder->SetOwner();\r
1653\r
1654 // add objects to folder\r
1655 while(i < size) {\r
1656 newFolder->Add(array->At(i));\r
1657 i++;\r
1658 }\r
1659 }\r
1660\r
1661return newFolder;\r
1662}\r
1663\r
1664//_____________________________________________________________________________\r
1665TFolder* AlidNdPtCorrection::CreateFolder(TString name,TString title) { \r
1666// create folder for analysed histograms\r
1667//\r
1668TFolder *folder = 0;\r
1669 folder = new TFolder(name.Data(),title.Data());\r
1670\r
1671 return folder;\r
1672}\r