]>
Commit | Line | Data |
---|---|---|
d6f334b5 | 1 | // $Id$ |
2 | // | |
3 | // Hadronic correction task. | |
4 | // | |
68263112 | 5 | // Author: R.Reed, C.Loizides |
d6f334b5 | 6 | |
7 | #include <TChain.h> | |
0b777a09 | 8 | #include <TClonesArray.h> |
35789a2d | 9 | #include <TH1F.h> |
10 | #include <TH2F.h> | |
d6f334b5 | 11 | #include <TList.h> |
35789a2d | 12 | #include <TLorentzVector.h> |
bcf52e4a | 13 | |
d6f334b5 | 14 | #include "AliAODCaloCluster.h" |
15 | #include "AliAODEvent.h" | |
0b777a09 | 16 | #include "AliAnalysisManager.h" |
d6f334b5 | 17 | #include "AliCentrality.h" |
18 | #include "AliESDCaloCluster.h" | |
0b777a09 | 19 | #include "AliESDtrack.h" |
0b777a09 | 20 | #include "AliPicoTrack.h" |
35789a2d | 21 | #include "AliVEventHandler.h" |
bcf52e4a | 22 | |
35789a2d | 23 | #include "AliHadCorrTask.h" |
0b777a09 | 24 | |
25 | ClassImp(AliHadCorrTask) | |
26 | ||
27 | //________________________________________________________________________ | |
d6f334b5 | 28 | AliHadCorrTask::AliHadCorrTask() : |
0b777a09 | 29 | AliAnalysisTaskSE("AliHadCorrTask"), |
d6f334b5 | 30 | fTracksName(), |
31 | fCaloName(), | |
32 | fOutCaloName(), | |
a3f43e7e | 33 | fPhiMatch(0.05), |
34 | fEtaMatch(0.025), | |
32bf39af | 35 | fDoTrackClus(0), |
f22cb876 | 36 | fHadCorr(0), |
0b777a09 | 37 | fMinPt(0.15), |
bcf52e4a | 38 | fCreateHisto(kFALSE), |
f22cb876 | 39 | fOutClusters(0), |
0b777a09 | 40 | fOutputList(0), |
0b777a09 | 41 | fHistNclusvsCent(0), |
42 | fHistNclusMatchvsCent(0), | |
f22cb876 | 43 | fHistEbefore(0), |
4711c521 | 44 | fHistEafter(0), |
45 | fHistEoPCent(0), | |
32bf39af | 46 | fHistNMatchCent(0), |
47 | fHistNMatchCent_trk(0), | |
48 | fHistCentrality(0) | |
0b777a09 | 49 | { |
d6f334b5 | 50 | // Default constructor. |
0b777a09 | 51 | |
a3f43e7e | 52 | for(Int_t i=0; i<8; i++) { |
53 | if (i<4) { | |
54 | fHistMatchEvsP[i] = 0; | |
55 | fHistMatchdRvsEP[i] = 0; | |
fbc9afc4 | 56 | for(Int_t j=0; j<3; j++) { |
57 | fHistEsubPch[i][j] = 0; | |
58 | fHistEsubPchRat[i][j] = 0; | |
59 | } | |
a3f43e7e | 60 | } |
fbc9afc4 | 61 | for(Int_t j=0; j<9; j++) { |
d6f334b5 | 62 | fHistMatchEtaPhi[i][j] = 0; |
a3f43e7e | 63 | } |
fbc9afc4 | 64 | |
a3f43e7e | 65 | } |
0b777a09 | 66 | } |
67 | ||
d6f334b5 | 68 | //________________________________________________________________________ |
69 | AliHadCorrTask::AliHadCorrTask(const char *name) : | |
70 | AliAnalysisTaskSE(name), | |
0b777a09 | 71 | fTracksName("Tracks"), |
72 | fCaloName("CaloClusters"), | |
d6f334b5 | 73 | fOutCaloName("CaloClustersCorr"), |
a3f43e7e | 74 | fPhiMatch(0.05), |
75 | fEtaMatch(0.025), | |
32bf39af | 76 | fDoTrackClus(0), |
f22cb876 | 77 | fHadCorr(0), |
f22cb876 | 78 | fMinPt(0.15), |
bcf52e4a | 79 | fCreateHisto(kFALSE), |
f22cb876 | 80 | fOutClusters(0), |
81 | fOutputList(0), | |
82 | fHistNclusvsCent(0), | |
83 | fHistNclusMatchvsCent(0), | |
84 | fHistEbefore(0), | |
4711c521 | 85 | fHistEafter(0), |
86 | fHistEoPCent(0), | |
32bf39af | 87 | fHistNMatchCent(0), |
88 | fHistNMatchCent_trk(0), | |
89 | fHistCentrality(0) | |
90 | ||
0b777a09 | 91 | { |
92 | // Standard constructor. | |
93 | ||
a3f43e7e | 94 | for(Int_t i=0; i<8; i++) { |
95 | if (i<4) { | |
fbc9afc4 | 96 | fHistMatchEvsP[i] = 0; |
97 | fHistMatchdRvsEP[i] = 0; | |
98 | for(Int_t j=0; j<3; j++) { | |
99 | fHistEsubPch[i][j] = 0; | |
100 | fHistEsubPchRat[i][j] = 0; | |
101 | } | |
a3f43e7e | 102 | } |
fbc9afc4 | 103 | for(Int_t j=0; j<9; j++) { |
d6f334b5 | 104 | fHistMatchEtaPhi[i][j] = 0; |
a3f43e7e | 105 | } |
106 | } | |
d6f334b5 | 107 | |
0b777a09 | 108 | fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex."; |
5d845887 | 109 | } |
110 | ||
111 | //________________________________________________________________________ | |
112 | AliHadCorrTask::AliHadCorrTask(const char *name, Bool_t histo) : | |
113 | AliAnalysisTaskSE(name), | |
114 | fTracksName("Tracks"), | |
115 | fCaloName("CaloClusters"), | |
116 | fOutCaloName("CaloClustersCorr"), | |
117 | fPhiMatch(0.05), | |
118 | fEtaMatch(0.025), | |
119 | fDoTrackClus(0), | |
120 | fHadCorr(0), | |
121 | fMinPt(0.15), | |
122 | fCreateHisto(histo), | |
123 | fOutClusters(0), | |
124 | fOutputList(0), | |
125 | fHistNclusvsCent(0), | |
126 | fHistNclusMatchvsCent(0), | |
127 | fHistEbefore(0), | |
128 | fHistEafter(0), | |
129 | fHistEoPCent(0), | |
130 | fHistNMatchCent(0), | |
131 | fHistNMatchCent_trk(0), | |
132 | fHistCentrality(0) | |
0b777a09 | 133 | |
5d845887 | 134 | { |
135 | // Standard constructor. | |
136 | ||
137 | for(Int_t i=0; i<8; i++) { | |
138 | if (i<4) { | |
139 | fHistMatchEvsP[i] = 0; | |
140 | fHistMatchdRvsEP[i] = 0; | |
141 | for(Int_t j=0; j<3; j++) { | |
142 | fHistEsubPch[i][j] = 0; | |
143 | fHistEsubPchRat[i][j] = 0; | |
144 | } | |
145 | } | |
146 | for(Int_t j=0; j<9; j++) { | |
147 | fHistMatchEtaPhi[i][j] = 0; | |
148 | } | |
149 | } | |
150 | ||
151 | fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex."; | |
152 | ||
153 | if (fCreateHisto) | |
154 | DefineOutput(1, TList::Class()); | |
d6f334b5 | 155 | } |
0b777a09 | 156 | |
157 | //________________________________________________________________________ | |
158 | AliHadCorrTask::~AliHadCorrTask() | |
159 | { | |
160 | // Destructor | |
161 | } | |
162 | ||
d6f334b5 | 163 | //________________________________________________________________________ |
164 | Int_t AliHadCorrTask::GetCentBin(Double_t cent) const | |
165 | { | |
166 | // Get centrality bin. | |
167 | ||
168 | Int_t centbin = -1; | |
169 | if (cent>=0 && cent<10) | |
170 | centbin=0; | |
171 | else if (cent>=10 && cent<30) | |
172 | centbin=1; | |
173 | else if (cent>=30 && cent<50) | |
174 | centbin=2; | |
175 | else if (cent>=50 && cent<=100) | |
176 | centbin=3; | |
177 | ||
178 | return centbin; | |
179 | } | |
180 | ||
181 | //________________________________________________________________________ | |
182 | Int_t AliHadCorrTask::GetMomBin(Double_t p) const | |
183 | { | |
184 | // Get momenum bin. | |
185 | ||
186 | Int_t pbin=-1; | |
187 | if (p>=0 && p<0.5) | |
188 | pbin=0; | |
fbc9afc4 | 189 | else if (p>=0.5 && p<1.0) |
d6f334b5 | 190 | pbin=1; |
fbc9afc4 | 191 | else if (p>=1.0 && p<1.5) |
d6f334b5 | 192 | pbin=2; |
fbc9afc4 | 193 | else if (p>=1.5 && p<2.) |
d6f334b5 | 194 | pbin=3; |
fbc9afc4 | 195 | else if (p>=2. && p<3.) |
d6f334b5 | 196 | pbin=4; |
fbc9afc4 | 197 | else if (p>=3. && p<4.) |
198 | pbin=5; | |
199 | else if (p>=4. && p<5.) | |
200 | pbin=6; | |
201 | else if (p>=5. && p<8.) | |
202 | pbin=7; | |
203 | else if (p>=8.) | |
204 | pbin=8; | |
d6f334b5 | 205 | |
206 | return pbin; | |
207 | } | |
208 | ||
f6e8a928 | 209 | //________________________________________________________________________ |
210 | Double_t AliHadCorrTask::GetEtaSigma(Int_t pbin) const | |
211 | { | |
212 | ||
213 | Double_t EtaSigma[9]={0.0097,0.0075,0.0059,0.0055,0.0053,0.005,0.005,0.045,0.042}; | |
214 | return 2.0*EtaSigma[pbin]; | |
215 | } | |
216 | //________________________________________________________________________ | |
217 | Double_t AliHadCorrTask::GetPhiMean(Int_t pbin, Int_t centbin) const | |
218 | { | |
219 | ||
220 | if (centbin==0){ | |
221 | Double_t PhiMean[9]={0.036, | |
222 | 0.021, | |
223 | 0.0121, | |
224 | 0.0084, | |
225 | 0.0060, | |
226 | 0.0041, | |
227 | 0.0031, | |
228 | 0.0022, | |
229 | 0.001}; | |
230 | return PhiMean[pbin]; | |
231 | }else if(centbin==1){ | |
232 | Double_t PhiMean[9]={0.036, | |
233 | 0.021, | |
234 | 0.0121, | |
235 | 0.0084, | |
236 | 0.0060, | |
237 | 0.0041, | |
238 | 0.0031, | |
239 | 0.0022, | |
240 | 0.001}; | |
241 | return PhiMean[pbin]; | |
242 | }else if(centbin==2){ | |
243 | Double_t PhiMean[9]={0.036, | |
244 | 0.021, | |
245 | 0.0121, | |
246 | 0.0084, | |
247 | 0.0060, | |
248 | 0.0041, | |
249 | 0.0031, | |
250 | 0.0022, | |
251 | 0.001}; | |
252 | return PhiMean[pbin]; | |
253 | }else if(centbin==3){ | |
254 | Double_t PhiMean[9]={0.036, | |
255 | 0.021, | |
256 | 0.0121, | |
257 | 0.0084, | |
258 | 0.0060, | |
259 | 0.0041, | |
260 | 0.0031, | |
261 | 0.0022, | |
262 | 0.001}; | |
263 | return PhiMean[pbin]; | |
264 | }else if(centbin==4){ | |
265 | Double_t PhiMean[9]={0.036, | |
266 | 0.021, | |
267 | 0.0127, | |
268 | 0.0089, | |
269 | 0.0068, | |
270 | 0.0049, | |
271 | 0.0038, | |
272 | 0.0028, | |
273 | 0.0018}; | |
274 | return PhiMean[pbin]*(-1.); | |
275 | }else if(centbin==5){ | |
276 | Double_t PhiMean[9]={0.036, | |
277 | 0.021, | |
278 | 0.0127, | |
279 | 0.0089, | |
280 | 0.0068, | |
281 | 0.0048, | |
282 | 0.0038, | |
283 | 0.0028, | |
284 | 0.0018}; | |
285 | return PhiMean[pbin]*(-1.); | |
286 | }else if(centbin==6){ | |
287 | Double_t PhiMean[9]={0.036, | |
288 | 0.021, | |
289 | 0.0127, | |
290 | 0.0089, | |
291 | 0.0068, | |
292 | 0.0045, | |
293 | 0.0035, | |
294 | 0.0028, | |
295 | 0.0018}; | |
296 | return PhiMean[pbin]*(-1.); | |
297 | }else if(centbin==7){ | |
298 | Double_t PhiMean[9]={0.036, | |
299 | 0.021, | |
300 | 0.0127, | |
301 | 0.0089, | |
302 | 0.0068, | |
303 | 0.0043, | |
304 | 0.0035, | |
305 | 0.0028, | |
306 | 0.0018}; | |
307 | return PhiMean[pbin]*(-1.); | |
308 | } | |
309 | ||
310 | return 0; | |
311 | ||
312 | } | |
313 | //________________________________________________________________________ | |
314 | Double_t AliHadCorrTask::GetPhiSigma(Int_t pbin, Int_t centbin) const | |
315 | { | |
316 | ||
317 | if (centbin==0){ | |
318 | Double_t PhiSigma[9]={0.0221, | |
319 | 0.0128, | |
320 | 0.0074, | |
321 | 0.0064, | |
322 | 0.0059, | |
323 | 0.0055, | |
324 | 0.0052, | |
325 | 0.0049, | |
326 | 0.0045}; | |
327 | return 2.*PhiSigma[pbin]; | |
328 | }else if(centbin==1){ | |
329 | Double_t PhiSigma[9]={0.0217, | |
330 | 0.0120, | |
331 | 0.0076, | |
332 | 0.0066, | |
333 | 0.0062, | |
334 | 0.0058, | |
335 | 0.0054, | |
336 | 0.0054, | |
337 | 0.0045}; | |
338 | return 2.*PhiSigma[pbin]; | |
339 | }else if(centbin==2){ | |
340 | Double_t PhiSigma[9]={0.0211, | |
341 | 0.0124, | |
342 | 0.0080, | |
343 | 0.0070, | |
344 | 0.0067, | |
345 | 0.0061, | |
346 | 0.0059, | |
347 | 0.0054, | |
348 | 0.0047}; | |
349 | return 2.*PhiSigma[pbin]; | |
350 | }else if(centbin==3){ | |
351 | Double_t PhiSigma[9]={0.0215, | |
352 | 0.0124, | |
353 | 0.0082, | |
354 | 0.0073, | |
355 | 0.0069, | |
356 | 0.0064, | |
357 | 0.0060, | |
358 | 0.0055, | |
359 | 0.0047}; | |
360 | return 2.*PhiSigma[pbin]; | |
361 | }else if(centbin==4){ | |
362 | Double_t PhiSigma[9]={0.0199, | |
363 | 0.0108, | |
364 | 0.0072, | |
365 | 0.0071, | |
366 | 0.0060, | |
367 | 0.0055, | |
368 | 0.0052, | |
369 | 0.0049, | |
370 | 0.0045}; | |
371 | return 2.*PhiSigma[pbin]; | |
372 | }else if(centbin==5){ | |
373 | Double_t PhiSigma[9]={0.0200, | |
374 | 0.0110, | |
375 | 0.0074, | |
376 | 0.0071, | |
377 | 0.0064, | |
378 | 0.0059, | |
379 | 0.0055, | |
380 | 0.0052, | |
381 | 0.0045}; | |
382 | return 2.*PhiSigma[pbin]; | |
383 | }else if(centbin==6){ | |
384 | Double_t PhiSigma[9]={0.0202, | |
385 | 0.0113, | |
386 | 0.0077, | |
387 | 0.0071, | |
388 | 0.0069, | |
389 | 0.0064, | |
390 | 0.0060, | |
391 | 0.0055, | |
392 | 0.0050}; | |
393 | return 2.*PhiSigma[pbin]; | |
394 | }else if(centbin==7){ | |
395 | Double_t PhiSigma[9]={0.0205, | |
396 | 0.0113, | |
397 | 0.0080, | |
398 | 0.0074, | |
399 | 0.0078, | |
400 | 0.0067, | |
401 | 0.0062, | |
402 | 0.0055, | |
403 | 0.0050}; | |
404 | return 2.*PhiSigma[pbin]; | |
405 | } | |
406 | ||
407 | return 0; | |
408 | ||
409 | } | |
410 | ||
0b777a09 | 411 | //________________________________________________________________________ |
412 | void AliHadCorrTask::UserCreateOutputObjects() | |
413 | { | |
d6f334b5 | 414 | // Create my user objects. |
0b777a09 | 415 | |
35789a2d | 416 | AliVEventHandler* handler = AliAnalysisManager::GetAnalysisManager()->GetInputEventHandler(); |
35789a2d | 417 | if (!handler) { |
418 | AliError("Input handler not available!"); | |
419 | return; | |
420 | } | |
421 | ||
422 | if (handler->InheritsFrom("AliESDInputHandler")) { | |
423 | fOutClusters = new TClonesArray("AliESDCaloCluster"); | |
d6f334b5 | 424 | } else if (handler->InheritsFrom("AliAODInputHandler")) { |
35789a2d | 425 | fOutClusters = new TClonesArray("AliAODCaloCluster"); |
d6f334b5 | 426 | } else { |
35789a2d | 427 | AliError("Input handler not recognized!"); |
428 | return; | |
429 | } | |
0b777a09 | 430 | fOutClusters->SetName(fOutCaloName); |
bcf52e4a | 431 | |
5d845887 | 432 | if (!fCreateHisto) |
433 | return; | |
434 | ||
0b777a09 | 435 | OpenFile(1); |
436 | fOutputList = new TList(); | |
437 | fOutputList->SetOwner(); | |
438 | ||
bcf52e4a | 439 | if (!fCreateHisto) { |
440 | PostData(1, fOutputList); | |
441 | return; | |
442 | } | |
443 | ||
a3f43e7e | 444 | for(Int_t icent=0; icent<8; ++icent) { |
fbc9afc4 | 445 | for(Int_t ipt=0; ipt<9; ++ipt) { |
d6f334b5 | 446 | TString name(Form("fHistMatchEtaPhi_%i_%i",icent,ipt)); |
447 | fHistMatchEtaPhi[icent][ipt] = new TH2F(name, name, 400, -0.2, 0.2, 1600, -0.8, 0.8); | |
4711c521 | 448 | fOutputList->Add(fHistMatchEtaPhi[icent][ipt]); |
449 | } | |
450 | ||
a3f43e7e | 451 | if(icent<4){ |
452 | TString name(Form("fHistMatchEvsP_%i",icent)); | |
453 | fHistMatchEvsP[icent] = new TH2F(name, name, 400, 0., 200., 1000, 0., 10.); | |
454 | fOutputList->Add(fHistMatchEvsP[icent]); | |
0b777a09 | 455 | |
a3f43e7e | 456 | name = Form("fHistMatchdRvsEP_%i",icent); |
457 | fHistMatchdRvsEP[icent] = new TH2F(name, name, 1000, 0., 1., 1000, 0., 10.); | |
458 | fOutputList->Add(fHistMatchdRvsEP[icent]); | |
459 | ||
fbc9afc4 | 460 | for(Int_t itrk=0; itrk<3; ++itrk){ |
461 | name = Form("fHistEsubPch_%i_%i",icent,itrk); | |
462 | fHistEsubPch[icent][itrk]=new TH1F(name, name, 400, 0., 100.); | |
463 | fHistEsubPch[icent][itrk]->Sumw2(); | |
464 | fOutputList->Add(fHistEsubPch[icent][itrk]); | |
465 | ||
466 | name = Form("fHistEsubPchRat_%i_%i",icent,itrk); | |
467 | fHistEsubPchRat[icent][itrk]=new TH2F(name, name, 400, 0., 200., 1000, 0., 10.); | |
468 | fOutputList->Add(fHistEsubPchRat[icent][itrk]); | |
469 | } | |
a3f43e7e | 470 | } |
d6f334b5 | 471 | } |
0b777a09 | 472 | |
32bf39af | 473 | fHistCentrality = new TH1F("fHistCentrality", "Centrality", 100, 0, 100); |
474 | ||
d6f334b5 | 475 | fHistNclusvsCent = new TH1F("Nclusvscent", "NclusVsCent", 100, 0, 100); |
476 | fHistNclusMatchvsCent = new TH1F("NclusMatchvscent", "NclusMatchVsCent", 100, 0, 100); | |
477 | fHistEbefore = new TH1F("Ebefore", "Ebefore", 100, 0, 100); | |
478 | fHistEafter = new TH1F("Eafter", "Eafter", 100, 0, 100); | |
479 | fHistEoPCent = new TH2F("EoPCent", "EoPCent", 100, 0, 100, 1000, 0, 10); | |
480 | fHistNMatchCent = new TH2F("NMatchesCent", "NMatchesCent", 100, 0, 100, 101, -0.5, 100.5); | |
32bf39af | 481 | fHistNMatchCent_trk = new TH2F("NMatchesCent_trk", "NMatchesCent_trk", 100, 0, 100, 101, -0.5, 100.5); |
0b777a09 | 482 | fOutputList->Add(fHistNclusMatchvsCent); |
483 | fOutputList->Add(fHistNclusvsCent); | |
484 | fOutputList->Add(fHistEbefore); | |
485 | fOutputList->Add(fHistEafter); | |
4711c521 | 486 | fOutputList->Add(fHistEoPCent); |
487 | fOutputList->Add(fHistNMatchCent); | |
32bf39af | 488 | fOutputList->Add(fHistNMatchCent_trk); |
489 | fOutputList->Add(fHistCentrality); | |
0b777a09 | 490 | |
491 | PostData(1, fOutputList); | |
492 | } | |
493 | ||
494 | //________________________________________________________________________ | |
495 | void AliHadCorrTask::UserExec(Option_t *) | |
496 | { | |
d6f334b5 | 497 | // Execute per event. |
0b777a09 | 498 | |
d6f334b5 | 499 | // post output in event if not yet present |
0b777a09 | 500 | if (!(InputEvent()->FindListObject(fOutCaloName))) |
501 | InputEvent()->AddObject(fOutClusters); | |
d6f334b5 | 502 | |
503 | // delete output | |
504 | fOutClusters->Delete(); | |
0b777a09 | 505 | |
d6f334b5 | 506 | // esd or aod mode |
507 | Bool_t esdMode = kTRUE; | |
508 | if (dynamic_cast<AliAODEvent*>(InputEvent())) | |
509 | esdMode = kFALSE; | |
0b777a09 | 510 | |
d6f334b5 | 511 | // optimization in case autobranch loading is off |
0b777a09 | 512 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); |
d6f334b5 | 513 | if (fCaloName == "CaloClusters") |
514 | am->LoadBranch("CaloClusters"); | |
515 | if (fTracksName == "Tracks") | |
516 | am->LoadBranch("Tracks"); | |
32bf39af | 517 | am->LoadBranch("Centrality"); |
518 | ||
519 | TList *l = InputEvent()->GetList(); | |
520 | ||
d6f334b5 | 521 | // get centrality |
c7b78119 | 522 | Double_t cent = -1; |
32bf39af | 523 | |
d6f334b5 | 524 | AliCentrality *centrality = InputEvent()->GetCentrality() ; |
32bf39af | 525 | |
d6f334b5 | 526 | if (centrality) |
527 | cent = centrality->GetCentralityPercentile("V0M"); | |
0b777a09 | 528 | else |
d6f334b5 | 529 | cent=99; // probably pp data |
32bf39af | 530 | |
d6f334b5 | 531 | if (cent<0) { |
20a2c788 | 532 | AliWarning(Form("Centrality negative: %f, assuming 99", cent)); |
533 | cent = 99; | |
d6f334b5 | 534 | } |
32bf39af | 535 | |
d6f334b5 | 536 | Int_t centbin = GetCentBin(cent); |
0b777a09 | 537 | |
bcf52e4a | 538 | if (fCreateHisto) |
539 | fHistCentrality->Fill(cent); | |
f6e8a928 | 540 | |
d6f334b5 | 541 | // get input collections |
542 | TClonesArray *tracks = 0; | |
543 | TClonesArray *clus = 0; | |
a3f43e7e | 544 | |
0b777a09 | 545 | tracks = dynamic_cast<TClonesArray*>(l->FindObject(fTracksName)); |
546 | if (!tracks) { | |
547 | AliError(Form("Pointer to tracks %s == 0", fTracksName.Data() )); | |
548 | return; | |
549 | } | |
0b777a09 | 550 | clus = dynamic_cast<TClonesArray*>(l->FindObject(fCaloName)); |
551 | if (!clus) { | |
552 | AliError(Form("Pointer to clus %s == 0", fCaloName.Data() )); | |
553 | return; | |
554 | } | |
555 | ||
d6f334b5 | 556 | // get primary vertex |
557 | Double_t vertex[3] = {0, 0, 0}; | |
558 | InputEvent()->GetPrimaryVertex()->GetXYZ(vertex); | |
559 | ||
560 | // loop over clusters and tracks | |
561 | const Int_t Nclus = clus->GetEntries(); | |
562 | const Int_t Ntrks = tracks->GetEntries(); | |
a3f43e7e | 563 | |
32bf39af | 564 | if (fDoTrackClus) { |
565 | for(Int_t t = 0; t<Ntrks; ++t) { | |
f6e8a928 | 566 | AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(t)); |
32bf39af | 567 | if (!track) |
568 | continue; | |
569 | Int_t Nmatches = 0; | |
570 | Double_t dEtaMin = 1e9; | |
571 | Double_t dPhiMin = 1e9; | |
572 | Int_t imin = -1; | |
573 | for(Int_t i=0; i < Nclus; ++i) { | |
f6e8a928 | 574 | AliVCluster *c = dynamic_cast<AliVCluster*>(clus->At(i)); |
32bf39af | 575 | if (!c) |
576 | continue; | |
577 | Double_t etadiff=999; | |
578 | Double_t phidiff=999; | |
579 | AliPicoTrack::GetEtaPhiDiff(track,c,phidiff,etadiff); | |
580 | Double_t dR = TMath::Sqrt(etadiff*etadiff+phidiff*phidiff); | |
581 | Double_t dRmin = TMath::Sqrt(dEtaMin*dEtaMin+dPhiMin*dPhiMin); | |
582 | if(dR > 25) | |
583 | continue; | |
584 | if(dR<dRmin){ | |
585 | dEtaMin = etadiff; | |
586 | dPhiMin = phidiff; | |
587 | imin = i; | |
588 | } | |
589 | if (TMath::Abs(phidiff)<0.05 && TMath::Abs(etadiff)<0.025) { // pp cuts!!! | |
590 | Nmatches++; | |
591 | } | |
592 | } | |
bcf52e4a | 593 | |
594 | if (fCreateHisto) | |
595 | fHistNMatchCent_trk->Fill(cent,Nmatches); | |
32bf39af | 596 | |
597 | track->SetEMCALcluster(imin); | |
598 | } | |
f6e8a928 | 599 | |
32bf39af | 600 | } |
601 | ||
d6f334b5 | 602 | for (Int_t iClus = 0, clusCount=0; iClus < Nclus; ++iClus) { |
f6e8a928 | 603 | AliVCluster *c = dynamic_cast<AliVCluster*>(clus->At(iClus)); |
d6f334b5 | 604 | if (!c) |
605 | continue; | |
606 | if (!c->IsEMCAL()) | |
607 | continue; | |
608 | ||
609 | // make primary particle out of cluster | |
610 | TLorentzVector nPart; | |
611 | c->GetMomentum(nPart, vertex); | |
612 | ||
613 | Double_t etclus = nPart.Pt(); | |
614 | if (etclus<fMinPt) | |
615 | continue; | |
616 | ||
617 | Double_t energyclus = nPart.P(); | |
618 | ||
619 | // reset cluster/track matching | |
620 | c->SetEmcCpvDistance(-1); | |
621 | c->SetTrackDistance(999,999); | |
622 | ||
623 | // run cluster-track matching | |
624 | Int_t imin = -1; | |
625 | Double_t dEtaMin = 1e9; | |
626 | Double_t dPhiMin = 1e9; | |
cce6a19c | 627 | Double_t dRmin = 1e9; |
d6f334b5 | 628 | Double_t totalTrkP = 0.0; // count total track momentum |
629 | Int_t Nmatches = 0; // count total number of matches | |
630 | for (Int_t t = 0; t<Ntrks; ++t) { | |
bcf52e4a | 631 | |
f6e8a928 | 632 | AliVTrack *track = dynamic_cast<AliVTrack*>(tracks->At(t)); |
bcf52e4a | 633 | |
d6f334b5 | 634 | if (!track) |
0b777a09 | 635 | continue; |
d6f334b5 | 636 | if (!track->IsEMCAL()) |
0b777a09 | 637 | continue; |
bcf52e4a | 638 | if (track->Pt() < fMinPt) |
d6f334b5 | 639 | continue; |
bcf52e4a | 640 | |
641 | Double_t etadiff = 999; | |
642 | Double_t phidiff = 999; | |
643 | AliPicoTrack::GetEtaPhiDiff(track, c, phidiff, etadiff); | |
644 | Double_t dR = TMath::Sqrt(etadiff * etadiff + phidiff * phidiff); | |
d6f334b5 | 645 | if(dR<dRmin){ |
646 | dEtaMin = etadiff; | |
647 | dPhiMin = phidiff; | |
648 | dRmin = dR; | |
649 | imin = t; | |
650 | } | |
f6e8a928 | 651 | |
652 | Double_t mom = track->P(); | |
653 | Int_t mombin = GetMomBin(mom); | |
654 | Int_t centbinch = centbin; | |
bcf52e4a | 655 | if (track->Charge() == -1 || track->Charge() == 255) |
f6e8a928 | 656 | centbinch += 4; |
657 | ||
bcf52e4a | 658 | if (fCreateHisto) { |
659 | if (fHadCorr > 1 && mombin > -1) { | |
660 | fHistMatchEtaPhi[centbinch][mombin]->Fill(etadiff, phidiff); | |
661 | fHistMatchdRvsEP[centbin]->Fill(dR, energyclus / mom); | |
662 | } | |
0b777a09 | 663 | } |
f6e8a928 | 664 | |
bcf52e4a | 665 | Double_t EtaCut = 0.0; |
666 | Double_t PhiCutlo = 0.0; | |
667 | Double_t PhiCuthi = 0.0; | |
668 | if(fPhiMatch > 0){ | |
669 | PhiCutlo = -fPhiMatch; | |
670 | PhiCuthi = fPhiMatch; | |
671 | } | |
672 | else { | |
673 | PhiCutlo = GetPhiMean(mombin,centbinch)-GetPhiSigma(mombin,centbin); | |
674 | PhiCuthi = GetPhiMean(mombin,centbinch)+GetPhiSigma(mombin,centbin); | |
f6e8a928 | 675 | } |
676 | ||
bcf52e4a | 677 | if(fEtaMatch > 0){ |
678 | EtaCut = fEtaMatch; | |
679 | } | |
680 | else { | |
681 | EtaCut = GetEtaSigma(mombin); | |
f6e8a928 | 682 | } |
683 | ||
bcf52e4a | 684 | if ((phidiff < PhiCuthi && phidiff > PhiCutlo) && TMath::Abs(etadiff) < EtaCut) { |
685 | if((fDoTrackClus && (track->GetEMCALcluster()) == iClus) || !fDoTrackClus){ | |
32bf39af | 686 | ++Nmatches; |
687 | totalTrkP += track->P(); | |
688 | } | |
d6f334b5 | 689 | } |
690 | } | |
691 | ||
692 | // store closest track | |
693 | c->SetEmcCpvDistance(imin); | |
694 | c->SetTrackDistance(dPhiMin, dEtaMin); | |
695 | ||
bcf52e4a | 696 | if(fCreateHisto) { |
697 | fHistNclusvsCent->Fill(cent); | |
698 | fHistEbefore->Fill(cent, energyclus); | |
699 | fHistNMatchCent->Fill(cent, Nmatches); | |
700 | ||
701 | if(Nmatches > 0) | |
702 | fHistNclusMatchvsCent->Fill(cent); | |
703 | } | |
32bf39af | 704 | |
d6f334b5 | 705 | // apply correction / subtraction |
bcf52e4a | 706 | if (fHadCorr > 0) { |
d6f334b5 | 707 | // to subtract only the closest track set fHadCor to a % |
708 | // to subtract all tracks within the cut set fHadCor to %+1 | |
bcf52e4a | 709 | if (fHadCorr > 1) { |
710 | if (totalTrkP > 0) { | |
711 | Double_t EoP = energyclus / totalTrkP; | |
a3f43e7e | 712 | Double_t Esub = (fHadCorr-1)*totalTrkP; |
bcf52e4a | 713 | if (Esub > energyclus) |
a3f43e7e | 714 | Esub = energyclus; |
fbc9afc4 | 715 | |
bcf52e4a | 716 | if(fCreateHisto) { |
717 | fHistEoPCent->Fill(cent, EoP); | |
718 | fHistMatchEvsP[centbin]->Fill(energyclus, EoP); | |
719 | ||
720 | if (Nmatches == 1) fHistEsubPchRat[centbin][0]->Fill(totalTrkP,Esub/totalTrkP); | |
721 | else if(Nmatches == 2) fHistEsubPchRat[centbin][1]->Fill(totalTrkP,Esub/totalTrkP); | |
722 | else fHistEsubPchRat[centbin][2]->Fill(totalTrkP,Esub/totalTrkP); | |
723 | ||
724 | if(Nmatches==1) fHistEsubPch[centbin][0]->Fill(totalTrkP,Esub); | |
725 | else if(Nmatches==2) fHistEsubPch[centbin][1]->Fill(totalTrkP,Esub); | |
726 | else fHistEsubPch[centbin][2]->Fill(totalTrkP,Esub); | |
727 | } | |
d6f334b5 | 728 | } |
bcf52e4a | 729 | energyclus -= (fHadCorr - 1) * totalTrkP; |
730 | } | |
731 | else if (imin >= 0) { | |
f6e8a928 | 732 | AliVTrack *t = dynamic_cast<AliVTrack*>(tracks->At(imin)); |
d6f334b5 | 733 | if (t) { |
734 | Double_t mom = t->P(); | |
735 | Int_t mombin = GetMomBin(mom); | |
a3f43e7e | 736 | Int_t centbinch = centbin; |
bcf52e4a | 737 | if (t->Charge() == -1 || t->Charge() == 255) |
a3f43e7e | 738 | centbinch += 4; |
bcf52e4a | 739 | |
740 | if (fCreateHisto) { | |
741 | fHistMatchEtaPhi[centbinch][mombin]->Fill(dEtaMin, dPhiMin); | |
742 | ||
743 | if (mom > 0) { | |
744 | fHistMatchEvsP[centbin]->Fill(energyclus, energyclus / mom); | |
745 | fHistEoPCent->Fill(cent, energyclus / mom); | |
746 | fHistMatchdRvsEP[centbin]->Fill(dRmin, energyclus / mom); | |
747 | } | |
748 | } | |
749 | ||
750 | Double_t EtaCut = 0.0; | |
751 | Double_t PhiCutlo = 0.0; | |
752 | Double_t PhiCuthi = 0.0; | |
753 | if(fPhiMatch > 0) { | |
754 | PhiCutlo = -fPhiMatch; | |
755 | PhiCuthi = fPhiMatch; | |
756 | } | |
757 | else { | |
758 | PhiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, centbin); | |
759 | PhiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, centbin); | |
f6e8a928 | 760 | } |
761 | ||
bcf52e4a | 762 | if(fEtaMatch > 0) { |
763 | EtaCut = fEtaMatch; | |
764 | } | |
765 | else { | |
766 | EtaCut = GetEtaSigma(mombin); | |
f6e8a928 | 767 | } |
768 | ||
bcf52e4a | 769 | if ((dPhiMin<PhiCuthi && dPhiMin>PhiCutlo) && TMath::Abs(dEtaMin) < EtaCut) { |
f6e8a928 | 770 | |
bcf52e4a | 771 | if((fDoTrackClus && (t->GetEMCALcluster()) == iClus) || !fDoTrackClus){ |
772 | energyclus -= fHadCorr * t->P(); | |
32bf39af | 773 | } |
d6f334b5 | 774 | } |
775 | } | |
0b777a09 | 776 | } |
777 | } | |
bcf52e4a | 778 | |
779 | if (energyclus < 0) | |
d6f334b5 | 780 | energyclus = 0; |
d6f334b5 | 781 | |
bcf52e4a | 782 | if (fCreateHisto) |
783 | fHistEafter->Fill(cent, energyclus); | |
784 | ||
785 | if (energyclus > 0) { // create corrected cluster | |
d6f334b5 | 786 | AliVCluster *oc; |
787 | if (esdMode) { | |
788 | oc = new ((*fOutClusters)[clusCount]) AliESDCaloCluster(*(dynamic_cast<AliESDCaloCluster*>(c))); | |
789 | } else { | |
790 | oc = new ((*fOutClusters)[clusCount]) AliAODCaloCluster(*(dynamic_cast<AliAODCaloCluster*>(c))); | |
791 | } | |
792 | oc->SetE(energyclus); | |
793 | ++clusCount; | |
794 | } | |
0b777a09 | 795 | } |
bcf52e4a | 796 | |
797 | if (fCreateHisto) | |
798 | PostData(1, fOutputList); | |
0b777a09 | 799 | } |
4711c521 | 800 | |
0b777a09 | 801 | //________________________________________________________________________ |
802 | void AliHadCorrTask::Terminate(Option_t *) | |
803 | { | |
d6f334b5 | 804 | // Nothing to be done. |
0b777a09 | 805 | } |
806 |