]>
Commit | Line | Data |
---|---|---|
d6f334b5 | 1 | // $Id$ |
2 | // | |
3 | // Hadronic correction task. | |
4 | // | |
1f6c8e39 | 5 | // Author: R.Reed, C.Loizides, S.Aiola |
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> |
bcf52e4a | 12 | |
3a385c49 | 13 | #include "AliAnalysisManager.h" |
d6f334b5 | 14 | #include "AliAODEvent.h" |
1f6c8e39 | 15 | #include "AliAODCaloCluster.h" |
d6f334b5 | 16 | #include "AliESDCaloCluster.h" |
1f6c8e39 | 17 | #include "AliVTrack.h" |
0b777a09 | 18 | #include "AliPicoTrack.h" |
35789a2d | 19 | #include "AliVEventHandler.h" |
1f6c8e39 | 20 | #include "AliEmcalParticle.h" |
1b4db73a | 21 | #include "AliEMCALGeometry.h" |
bcf52e4a | 22 | |
35789a2d | 23 | #include "AliHadCorrTask.h" |
0b777a09 | 24 | |
25 | ClassImp(AliHadCorrTask) | |
26 | ||
27 | //________________________________________________________________________ | |
d6f334b5 | 28 | AliHadCorrTask::AliHadCorrTask() : |
1f6c8e39 | 29 | AliAnalysisTaskEmcal("AliHadCorrTask", kFALSE), |
d6f334b5 | 30 | fOutCaloName(), |
a3f43e7e | 31 | fPhiMatch(0.05), |
32 | fEtaMatch(0.025), | |
32bf39af | 33 | fDoTrackClus(0), |
f22cb876 | 34 | fHadCorr(0), |
1f6c8e39 | 35 | fEexclCell(0), |
8c0d179d | 36 | fDoExact(kFALSE), |
1b4db73a | 37 | fEsdMode(kTRUE), |
f22cb876 | 38 | fOutClusters(0), |
0b777a09 | 39 | fHistNclusvsCent(0), |
40 | fHistNclusMatchvsCent(0), | |
f22cb876 | 41 | fHistEbefore(0), |
4711c521 | 42 | fHistEafter(0), |
43 | fHistEoPCent(0), | |
32bf39af | 44 | fHistNMatchCent(0), |
1b4db73a | 45 | fHistNClusMatchCent(0) |
0b777a09 | 46 | { |
d6f334b5 | 47 | // Default constructor. |
0b777a09 | 48 | |
a3f43e7e | 49 | for(Int_t i=0; i<8; i++) { |
9f52c61f | 50 | fHistEsubPch[i] = 0; |
51 | fHistEsubPchRat[i] = 0; | |
0c5f6f08 | 52 | fHistEsubPchRatAll[i] = 0; |
9f52c61f | 53 | |
a3f43e7e | 54 | if (i<4) { |
1f6c8e39 | 55 | fHistMatchEvsP[i] = 0; |
56 | fHistMatchdRvsEP[i] = 0; | |
57 | fHistNMatchEnergy[i] = 0; | |
9f52c61f | 58 | |
8c0d179d | 59 | fHistEmbTrackMatchesOversub[i] = 0; |
60 | fHistNonEmbTrackMatchesOversub[i] = 0; | |
88646d1d | 61 | fHistOversubMCClusters[i] = 0; |
62 | fHistOversubNonMCClusters[i] = 0; | |
8c0d179d | 63 | fHistOversub[i] = 0; |
8c0d179d | 64 | |
9f52c61f | 65 | for(Int_t j=0; j<4; j++) |
66 | fHistNCellsEnergy[i][j] = 0; | |
a3f43e7e | 67 | } |
9f52c61f | 68 | |
fbc9afc4 | 69 | for(Int_t j=0; j<9; j++) { |
9f52c61f | 70 | for(Int_t k=0; k<2; k++) |
1f6c8e39 | 71 | fHistMatchEtaPhi[i][j][k] = 0; |
a3f43e7e | 72 | } |
73 | } | |
0b777a09 | 74 | } |
75 | ||
5d845887 | 76 | //________________________________________________________________________ |
77 | AliHadCorrTask::AliHadCorrTask(const char *name, Bool_t histo) : | |
1f6c8e39 | 78 | AliAnalysisTaskEmcal(name, histo), |
5d845887 | 79 | fOutCaloName("CaloClustersCorr"), |
80 | fPhiMatch(0.05), | |
81 | fEtaMatch(0.025), | |
a89e2005 | 82 | fDoTrackClus(1), |
5d845887 | 83 | fHadCorr(0), |
1f6c8e39 | 84 | fEexclCell(0), |
8c0d179d | 85 | fDoExact(kFALSE), |
1b4db73a | 86 | fEsdMode(kTRUE), |
5d845887 | 87 | fOutClusters(0), |
5d845887 | 88 | fHistNclusvsCent(0), |
89 | fHistNclusMatchvsCent(0), | |
90 | fHistEbefore(0), | |
91 | fHistEafter(0), | |
92 | fHistEoPCent(0), | |
93 | fHistNMatchCent(0), | |
1b4db73a | 94 | fHistNClusMatchCent(0) |
5d845887 | 95 | { |
96 | // Standard constructor. | |
97 | ||
1b4db73a | 98 | for(Int_t i=0; i<8; i++) { |
9f52c61f | 99 | fHistEsubPch[i] = 0; |
100 | fHistEsubPchRat[i] = 0; | |
0c5f6f08 | 101 | fHistEsubPchRatAll[i] = 0; |
9f52c61f | 102 | |
5d845887 | 103 | if (i<4) { |
9f52c61f | 104 | fHistMatchEvsP[i] = 0; |
105 | fHistMatchdRvsEP[i] = 0; | |
106 | fHistNMatchEnergy[i] = 0; | |
107 | ||
8c0d179d | 108 | fHistEmbTrackMatchesOversub[i] = 0; |
109 | fHistNonEmbTrackMatchesOversub[i] = 0; | |
88646d1d | 110 | fHistOversubMCClusters[i] = 0; |
111 | fHistOversubNonMCClusters[i] = 0; | |
8c0d179d | 112 | fHistOversub[i] = 0; |
8c0d179d | 113 | |
9f52c61f | 114 | for(Int_t j=0; j<4; j++) |
1f6c8e39 | 115 | fHistNCellsEnergy[i][j] = 0; |
5d845887 | 116 | } |
9f52c61f | 117 | |
5d845887 | 118 | for(Int_t j=0; j<9; j++) { |
9f52c61f | 119 | for(Int_t k=0; k<2; k++) |
1f6c8e39 | 120 | fHistMatchEtaPhi[i][j][k] = 0; |
5d845887 | 121 | } |
122 | } | |
1b4db73a | 123 | |
124 | SetMakeGeneralHistograms(histo); | |
5d845887 | 125 | |
126 | fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex."; | |
d6f334b5 | 127 | } |
0b777a09 | 128 | |
129 | //________________________________________________________________________ | |
130 | AliHadCorrTask::~AliHadCorrTask() | |
131 | { | |
132 | // Destructor | |
133 | } | |
134 | ||
d6f334b5 | 135 | //________________________________________________________________________ |
b8b1b9e1 | 136 | UInt_t AliHadCorrTask::GetMomBin(Double_t p) const |
d6f334b5 | 137 | { |
138 | // Get momenum bin. | |
139 | ||
b8b1b9e1 | 140 | UInt_t pbin=0; |
4f833c06 | 141 | if (p<0.5) |
d6f334b5 | 142 | pbin=0; |
fbc9afc4 | 143 | else if (p>=0.5 && p<1.0) |
d6f334b5 | 144 | pbin=1; |
fbc9afc4 | 145 | else if (p>=1.0 && p<1.5) |
d6f334b5 | 146 | pbin=2; |
fbc9afc4 | 147 | else if (p>=1.5 && p<2.) |
d6f334b5 | 148 | pbin=3; |
fbc9afc4 | 149 | else if (p>=2. && p<3.) |
d6f334b5 | 150 | pbin=4; |
fbc9afc4 | 151 | else if (p>=3. && p<4.) |
152 | pbin=5; | |
153 | else if (p>=4. && p<5.) | |
154 | pbin=6; | |
155 | else if (p>=5. && p<8.) | |
156 | pbin=7; | |
b9a9ae7e | 157 | else |
fbc9afc4 | 158 | pbin=8; |
d6f334b5 | 159 | |
160 | return pbin; | |
161 | } | |
162 | ||
f6e8a928 | 163 | //________________________________________________________________________ |
164 | Double_t AliHadCorrTask::GetEtaSigma(Int_t pbin) const | |
165 | { | |
4f833c06 | 166 | // Get sigma in eta. |
167 | ||
1f6c8e39 | 168 | Double_t EtaSigma[9]={0.0097,0.0075,0.0059,0.0055,0.0053,0.005,0.005,0.0045,0.0042}; |
f6e8a928 | 169 | return 2.0*EtaSigma[pbin]; |
170 | } | |
1f6c8e39 | 171 | |
f6e8a928 | 172 | //________________________________________________________________________ |
173 | Double_t AliHadCorrTask::GetPhiMean(Int_t pbin, Int_t centbin) const | |
174 | { | |
4f833c06 | 175 | // Get phi mean. |
176 | ||
177 | if (centbin==0) { | |
f6e8a928 | 178 | Double_t PhiMean[9]={0.036, |
179 | 0.021, | |
180 | 0.0121, | |
181 | 0.0084, | |
182 | 0.0060, | |
183 | 0.0041, | |
184 | 0.0031, | |
185 | 0.0022, | |
186 | 0.001}; | |
187 | return PhiMean[pbin]; | |
4f833c06 | 188 | } else if (centbin==1) { |
f6e8a928 | 189 | Double_t PhiMean[9]={0.036, |
190 | 0.021, | |
191 | 0.0121, | |
192 | 0.0084, | |
193 | 0.0060, | |
194 | 0.0041, | |
195 | 0.0031, | |
196 | 0.0022, | |
197 | 0.001}; | |
198 | return PhiMean[pbin]; | |
4f833c06 | 199 | } else if (centbin==2) { |
f6e8a928 | 200 | Double_t PhiMean[9]={0.036, |
201 | 0.021, | |
202 | 0.0121, | |
203 | 0.0084, | |
204 | 0.0060, | |
205 | 0.0041, | |
206 | 0.0031, | |
207 | 0.0022, | |
208 | 0.001}; | |
209 | return PhiMean[pbin]; | |
4f833c06 | 210 | } else if (centbin==3) { |
f6e8a928 | 211 | Double_t PhiMean[9]={0.036, |
212 | 0.021, | |
213 | 0.0121, | |
214 | 0.0084, | |
215 | 0.0060, | |
216 | 0.0041, | |
217 | 0.0031, | |
218 | 0.0022, | |
219 | 0.001}; | |
220 | return PhiMean[pbin]; | |
4f833c06 | 221 | } else if (centbin==4) { |
f6e8a928 | 222 | Double_t PhiMean[9]={0.036, |
223 | 0.021, | |
224 | 0.0127, | |
225 | 0.0089, | |
226 | 0.0068, | |
227 | 0.0049, | |
228 | 0.0038, | |
229 | 0.0028, | |
230 | 0.0018}; | |
231 | return PhiMean[pbin]*(-1.); | |
4f833c06 | 232 | } else if (centbin==5) { |
f6e8a928 | 233 | Double_t PhiMean[9]={0.036, |
234 | 0.021, | |
235 | 0.0127, | |
236 | 0.0089, | |
237 | 0.0068, | |
238 | 0.0048, | |
239 | 0.0038, | |
240 | 0.0028, | |
241 | 0.0018}; | |
242 | return PhiMean[pbin]*(-1.); | |
4f833c06 | 243 | } else if (centbin==6) { |
f6e8a928 | 244 | Double_t PhiMean[9]={0.036, |
245 | 0.021, | |
246 | 0.0127, | |
247 | 0.0089, | |
248 | 0.0068, | |
249 | 0.0045, | |
250 | 0.0035, | |
251 | 0.0028, | |
252 | 0.0018}; | |
253 | return PhiMean[pbin]*(-1.); | |
4f833c06 | 254 | } else if (centbin==7) { |
f6e8a928 | 255 | Double_t PhiMean[9]={0.036, |
256 | 0.021, | |
257 | 0.0127, | |
258 | 0.0089, | |
259 | 0.0068, | |
260 | 0.0043, | |
261 | 0.0035, | |
262 | 0.0028, | |
263 | 0.0018}; | |
264 | return PhiMean[pbin]*(-1.); | |
265 | } | |
266 | ||
267 | return 0; | |
f6e8a928 | 268 | } |
1f6c8e39 | 269 | |
f6e8a928 | 270 | //________________________________________________________________________ |
271 | Double_t AliHadCorrTask::GetPhiSigma(Int_t pbin, Int_t centbin) const | |
272 | { | |
4f833c06 | 273 | // Get phi sigma. |
274 | ||
275 | if (centbin==0) { | |
f6e8a928 | 276 | Double_t PhiSigma[9]={0.0221, |
277 | 0.0128, | |
278 | 0.0074, | |
279 | 0.0064, | |
280 | 0.0059, | |
281 | 0.0055, | |
282 | 0.0052, | |
283 | 0.0049, | |
284 | 0.0045}; | |
285 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 286 | } else if (centbin==1) { |
f6e8a928 | 287 | Double_t PhiSigma[9]={0.0217, |
288 | 0.0120, | |
289 | 0.0076, | |
290 | 0.0066, | |
291 | 0.0062, | |
292 | 0.0058, | |
293 | 0.0054, | |
294 | 0.0054, | |
295 | 0.0045}; | |
296 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 297 | } else if (centbin==2) { |
f6e8a928 | 298 | Double_t PhiSigma[9]={0.0211, |
299 | 0.0124, | |
300 | 0.0080, | |
301 | 0.0070, | |
302 | 0.0067, | |
303 | 0.0061, | |
304 | 0.0059, | |
305 | 0.0054, | |
306 | 0.0047}; | |
307 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 308 | } else if (centbin==3) { |
f6e8a928 | 309 | Double_t PhiSigma[9]={0.0215, |
310 | 0.0124, | |
311 | 0.0082, | |
312 | 0.0073, | |
313 | 0.0069, | |
314 | 0.0064, | |
315 | 0.0060, | |
316 | 0.0055, | |
317 | 0.0047}; | |
318 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 319 | } else if (centbin==4) { |
f6e8a928 | 320 | Double_t PhiSigma[9]={0.0199, |
321 | 0.0108, | |
322 | 0.0072, | |
323 | 0.0071, | |
324 | 0.0060, | |
325 | 0.0055, | |
326 | 0.0052, | |
327 | 0.0049, | |
328 | 0.0045}; | |
329 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 330 | } else if (centbin==5) { |
f6e8a928 | 331 | Double_t PhiSigma[9]={0.0200, |
332 | 0.0110, | |
333 | 0.0074, | |
334 | 0.0071, | |
335 | 0.0064, | |
336 | 0.0059, | |
337 | 0.0055, | |
338 | 0.0052, | |
339 | 0.0045}; | |
340 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 341 | } else if (centbin==6) { |
f6e8a928 | 342 | Double_t PhiSigma[9]={0.0202, |
343 | 0.0113, | |
344 | 0.0077, | |
345 | 0.0071, | |
346 | 0.0069, | |
347 | 0.0064, | |
348 | 0.0060, | |
349 | 0.0055, | |
350 | 0.0050}; | |
351 | return 2.*PhiSigma[pbin]; | |
4f833c06 | 352 | } else if (centbin==7) { |
f6e8a928 | 353 | Double_t PhiSigma[9]={0.0205, |
354 | 0.0113, | |
355 | 0.0080, | |
356 | 0.0074, | |
357 | 0.0078, | |
358 | 0.0067, | |
359 | 0.0062, | |
360 | 0.0055, | |
361 | 0.0050}; | |
362 | return 2.*PhiSigma[pbin]; | |
363 | } | |
364 | ||
365 | return 0; | |
f6e8a928 | 366 | } |
367 | ||
0b777a09 | 368 | //________________________________________________________________________ |
369 | void AliHadCorrTask::UserCreateOutputObjects() | |
370 | { | |
d6f334b5 | 371 | // Create my user objects. |
0b777a09 | 372 | |
5d845887 | 373 | if (!fCreateHisto) |
374 | return; | |
375 | ||
1b4db73a | 376 | AliAnalysisTaskEmcal::UserCreateOutputObjects(); |
0b777a09 | 377 | |
1f6c8e39 | 378 | TString name; |
bcf52e4a | 379 | |
8c0d179d | 380 | const Int_t nCentChBins = fNcentBins * 2; |
381 | ||
382 | for(Int_t icent=0; icent<nCentChBins; ++icent) { | |
fbc9afc4 | 383 | for(Int_t ipt=0; ipt<9; ++ipt) { |
1f6c8e39 | 384 | for(Int_t ieta=0; ieta<2; ++ieta) { |
1f6c8e39 | 385 | name = Form("fHistMatchEtaPhi_%i_%i_%i",icent,ipt,ieta); |
8c0d179d | 386 | fHistMatchEtaPhi[icent][ipt][ieta] = new TH2F(name, name, fNbins, -0.1, 0.1, fNbins, -0.1, 0.1); |
1f6c8e39 | 387 | fOutput->Add(fHistMatchEtaPhi[icent][ipt][ieta]); |
388 | } | |
4711c521 | 389 | } |
390 | ||
1f6c8e39 | 391 | name = Form("fHistEsubPch_%i",icent); |
8c0d179d | 392 | fHistEsubPch[icent]=new TH1F(name, name, fNbins, fMinBinPt, fMaxBinPt); |
1f6c8e39 | 393 | fOutput->Add(fHistEsubPch[icent]); |
394 | ||
395 | name = Form("fHistEsubPchRat_%i",icent); | |
8c0d179d | 396 | fHistEsubPchRat[icent]=new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.); |
1f6c8e39 | 397 | fOutput->Add(fHistEsubPchRat[icent]); |
0c5f6f08 | 398 | |
399 | name = Form("fHistEsubPchRatAll_%i",icent); | |
8c0d179d | 400 | fHistEsubPchRatAll[icent]=new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.); |
0c5f6f08 | 401 | fOutput->Add(fHistEsubPchRatAll[icent]); |
1f6c8e39 | 402 | |
8c0d179d | 403 | if (icent<fNcentBins) { |
1b4db73a | 404 | for(Int_t itrk=0; itrk<4; ++itrk) { |
405 | name = Form("fHistNCellsEnergy_%i_%i",icent,itrk); | |
8c0d179d | 406 | fHistNCellsEnergy[icent][itrk] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, 101, -0.5, 100.5); |
1b4db73a | 407 | fOutput->Add(fHistNCellsEnergy[icent][itrk]); |
408 | } | |
409 | ||
1f6c8e39 | 410 | name = Form("fHistMatchEvsP_%i",icent); |
8c0d179d | 411 | fHistMatchEvsP[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.); |
1f6c8e39 | 412 | fOutput->Add(fHistMatchEvsP[icent]); |
0b777a09 | 413 | |
a3f43e7e | 414 | name = Form("fHistMatchdRvsEP_%i",icent); |
8c0d179d | 415 | fHistMatchdRvsEP[icent] = new TH2F(name, name, fNbins, 0., 0.2, fNbins*2, 0., 10.); |
1f6c8e39 | 416 | fOutput->Add(fHistMatchdRvsEP[icent]); |
417 | ||
418 | name = Form("fHistNMatchEnergy_%i",icent); | |
8c0d179d | 419 | fHistNMatchEnergy[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, 101, -0.5, 100.5); |
1f6c8e39 | 420 | fOutput->Add(fHistNMatchEnergy[icent]); |
a3f43e7e | 421 | } |
d6f334b5 | 422 | } |
0b777a09 | 423 | |
d6f334b5 | 424 | fHistNclusvsCent = new TH1F("Nclusvscent", "NclusVsCent", 100, 0, 100); |
425 | fHistNclusMatchvsCent = new TH1F("NclusMatchvscent", "NclusMatchVsCent", 100, 0, 100); | |
8c0d179d | 426 | fHistEbefore = new TH1F("Ebefore", "Ebefore", fNbins, fMinBinPt, fMaxBinPt); |
427 | fHistEafter = new TH1F("Eafter", "Eafter", fNbins, fMinBinPt, fMaxBinPt); | |
428 | fHistEoPCent = new TH2F("EoPCent", "EoPCent", 100, 0, 100, fNbins*2, 0, 10); | |
429 | fHistNMatchCent = new TH2F("NMatchesCent", "NMatchesCent", 100, 0, 100, 11, -0.5, 10.5); | |
430 | fHistNClusMatchCent = new TH2F("NClusMatchesCent", "NClusMatchesCent", 100, 0, 100, 11, -0.5, 10.5); | |
1f6c8e39 | 431 | |
432 | fOutput->Add(fHistNclusMatchvsCent); | |
433 | fOutput->Add(fHistNclusvsCent); | |
434 | fOutput->Add(fHistEbefore); | |
435 | fOutput->Add(fHistEafter); | |
436 | fOutput->Add(fHistEoPCent); | |
437 | fOutput->Add(fHistNMatchCent); | |
1b4db73a | 438 | fOutput->Add(fHistNClusMatchCent); |
1f6c8e39 | 439 | |
8c0d179d | 440 | if (fIsEmbedded) { |
441 | for(Int_t icent=0; icent<fNcentBins; ++icent) { | |
442 | name = Form("fHistEmbTrackMatchesOversub_%d",icent); | |
88646d1d | 443 | fHistEmbTrackMatchesOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2); |
8c0d179d | 444 | fHistEmbTrackMatchesOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)"); |
88646d1d | 445 | fHistEmbTrackMatchesOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}"); |
8c0d179d | 446 | fOutput->Add(fHistEmbTrackMatchesOversub[icent]); |
447 | ||
448 | name = Form("fHistNonEmbTrackMatchesOversub_%d",icent); | |
88646d1d | 449 | fHistNonEmbTrackMatchesOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2); |
8c0d179d | 450 | fHistNonEmbTrackMatchesOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)"); |
88646d1d | 451 | fHistNonEmbTrackMatchesOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}"); |
8c0d179d | 452 | fOutput->Add(fHistNonEmbTrackMatchesOversub[icent]); |
453 | ||
88646d1d | 454 | name = Form("fHistOversubMCClusters_%d",icent); |
455 | fHistOversubMCClusters[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2); | |
456 | fHistOversubMCClusters[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)"); | |
457 | fHistOversubMCClusters[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}"); | |
458 | fOutput->Add(fHistOversubMCClusters[icent]); | |
459 | ||
460 | name = Form("fHistOversubNonMCClusters_%d",icent); | |
461 | fHistOversubNonMCClusters[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2); | |
462 | fHistOversubNonMCClusters[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)"); | |
463 | fHistOversubNonMCClusters[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}"); | |
464 | fOutput->Add(fHistOversubNonMCClusters[icent]); | |
465 | ||
8c0d179d | 466 | name = Form("fHistOversub_%d",icent); |
88646d1d | 467 | fHistOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2); |
8c0d179d | 468 | fHistOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)"); |
88646d1d | 469 | fHistOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}"); |
8c0d179d | 470 | fOutput->Add(fHistOversub[icent]); |
8c0d179d | 471 | } |
472 | } | |
473 | ||
1f6c8e39 | 474 | PostData(1, fOutput); |
0b777a09 | 475 | } |
476 | ||
1b4db73a | 477 | //________________________________________________________________________ |
478 | void AliHadCorrTask::ExecOnce() | |
479 | { | |
0983dd5b | 480 | // Do base class initializations and if it fails -> bail out |
481 | if (!fInitialized) | |
482 | AliAnalysisTaskEmcal::ExecOnce(); | |
7030f36f | 483 | |
484 | if (!fInitialized) | |
485 | return; | |
486 | ||
1b4db73a | 487 | if (dynamic_cast<AliAODEvent*>(InputEvent())) |
488 | fEsdMode = kFALSE; | |
489 | ||
490 | if (fEsdMode) { // optimization in case autobranch loading is off | |
491 | AliAnalysisManager *am = AliAnalysisManager::GetAnalysisManager(); | |
492 | if (fCaloName == "CaloClusters") | |
493 | am->LoadBranch("CaloClusters"); | |
494 | if (fTracksName == "Tracks") | |
495 | am->LoadBranch("Tracks"); | |
496 | am->LoadBranch("Centrality."); | |
497 | } | |
498 | ||
499 | if (fEsdMode) | |
500 | fOutClusters = new TClonesArray("AliESDCaloCluster"); | |
501 | else | |
502 | fOutClusters = new TClonesArray("AliAODCaloCluster"); | |
503 | ||
504 | fOutClusters->SetName(fOutCaloName); | |
505 | ||
506 | // post output in event if not yet present | |
507 | if (!(InputEvent()->FindListObject(fOutCaloName))) { | |
508 | InputEvent()->AddObject(fOutClusters); | |
509 | } | |
510 | else { | |
7030f36f | 511 | fInitialized = kFALSE; |
1b4db73a | 512 | AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutCaloName.Data())); |
513 | return; | |
514 | } | |
1b4db73a | 515 | } |
516 | ||
517 | //________________________________________________________________________ | |
518 | void AliHadCorrTask::DoTrackLoop() | |
519 | { | |
520 | const Int_t Ntracks = fTracks->GetEntries(); | |
521 | for (Int_t iTrk = 0; iTrk < Ntracks; ++iTrk) { | |
522 | Int_t NmatchClus = 0; | |
523 | ||
524 | AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iTrk)); | |
525 | if (!emctrack) | |
526 | continue; | |
527 | if (!emctrack->IsEMCAL()) | |
528 | continue; | |
529 | ||
530 | AliVTrack *track = emctrack->GetTrack(); | |
531 | if (!track) | |
532 | continue; | |
533 | if (!AcceptTrack(track)) | |
534 | continue; | |
535 | ||
43032ce2 | 536 | if (track->GetTrackEtaOnEMCal() < fGeom->GetArm1EtaMin() - fEtaMatch || |
537 | track->GetTrackEtaOnEMCal() > fGeom->GetArm1EtaMax() + fEtaMatch || | |
538 | track->GetTrackPhiOnEMCal() < fGeom->GetArm1PhiMin() * TMath::DegToRad() - fPhiMatch || | |
539 | track->GetTrackPhiOnEMCal() > fGeom->GetArm1PhiMax() * TMath::DegToRad() + fPhiMatch) | |
1b4db73a | 540 | continue; |
541 | ||
542 | Int_t Nclus = emctrack->GetNumberOfMatchedObj(); | |
543 | ||
544 | for (Int_t iClus = 0; iClus < Nclus; ++iClus) { | |
545 | AliEmcalParticle *emccluster = static_cast<AliEmcalParticle*>(fCaloClusters->At(iClus)); | |
546 | if (!emccluster) | |
547 | continue; | |
548 | ||
549 | AliVCluster *cluster = emccluster->GetCluster(); | |
550 | if (!cluster) | |
551 | continue; | |
552 | if (!AcceptCluster(cluster)) | |
553 | continue; | |
554 | ||
555 | Double_t etadiff = 999; | |
556 | Double_t phidiff = 999; | |
557 | AliPicoTrack::GetEtaPhiDiff(track, cluster, phidiff, etadiff); | |
558 | ||
559 | if (TMath::Abs(phidiff) < fPhiMatch && TMath::Abs(etadiff) < fEtaMatch) | |
560 | NmatchClus++; | |
561 | } | |
562 | fHistNClusMatchCent->Fill(fCent, NmatchClus); | |
563 | } | |
564 | } | |
565 | ||
1f6c8e39 | 566 | //________________________________________________________________________ |
8c0d179d | 567 | void AliHadCorrTask::DoMatchedTracksLoop(AliEmcalParticle *emccluster, Double_t &totalTrkP, Int_t &Nmatches, Double_t &trkPMCfrac, Int_t &NMCmatches) |
1f6c8e39 | 568 | { |
c7d52659 | 569 | // Do the loop over matched tracks for cluster emccluster. |
1f6c8e39 | 570 | |
571 | AliVCluster *cluster = emccluster->GetCluster(); | |
572 | Int_t iClus = emccluster->IdInCollection(); | |
573 | Double_t energyclus = cluster->E(); | |
574 | ||
575 | // loop over matched tracks | |
c7d52659 | 576 | const Int_t Ntrks = emccluster->GetNumberOfMatchedObj(); |
1f6c8e39 | 577 | for (Int_t i = 0; i < Ntrks; ++i) { |
578 | Int_t iTrack = emccluster->GetMatchedObjId(i); | |
579 | Double_t dR = emccluster->GetMatchedObjDistance(i); | |
580 | ||
4f833c06 | 581 | AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iTrack)); |
1f6c8e39 | 582 | if (!emctrack) |
583 | continue; | |
56bd3193 | 584 | |
585 | // check if track also points to cluster | |
586 | if (fDoTrackClus && (emctrack->GetMatchedObjId(0)) != iClus) | |
587 | continue; | |
588 | ||
1f6c8e39 | 589 | AliVTrack *track = emctrack->GetTrack(); |
590 | if (!track) | |
591 | continue; | |
592 | if (!AcceptTrack(track)) | |
593 | continue; | |
594 | ||
595 | Double_t etadiff = 999; | |
596 | Double_t phidiff = 999; | |
597 | AliPicoTrack::GetEtaPhiDiff(track, cluster, phidiff, etadiff); | |
56bd3193 | 598 | |
1f6c8e39 | 599 | Double_t mom = track->P(); |
b8b1b9e1 | 600 | UInt_t mombin = GetMomBin(mom); |
1f6c8e39 | 601 | Int_t centbinch = fCentBin; |
c7d52659 | 602 | if (track->Charge()<0) |
8c0d179d | 603 | centbinch += fNcentBins; |
1f6c8e39 | 604 | |
605 | if (fCreateHisto) { | |
606 | Int_t etabin = 0; | |
607 | if(track->Eta() > 0) etabin=1; | |
1f6c8e39 | 608 | fHistMatchEtaPhi[centbinch][mombin][etabin]->Fill(etadiff, phidiff); |
609 | } | |
610 | ||
c7d52659 | 611 | Double_t etaCut = 0.0; |
612 | Double_t phiCutlo = 0.0; | |
613 | Double_t phiCuthi = 0.0; | |
1f6c8e39 | 614 | |
4f833c06 | 615 | if (fPhiMatch > 0) { |
c7d52659 | 616 | phiCutlo = -fPhiMatch; |
875b35fc | 617 | phiCuthi = +fPhiMatch; |
618 | } else { | |
c7d52659 | 619 | phiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, fCentBin); |
620 | phiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, fCentBin); | |
1f6c8e39 | 621 | } |
622 | ||
623 | if (fEtaMatch > 0) { | |
c7d52659 | 624 | etaCut = fEtaMatch; |
875b35fc | 625 | } else { |
c7d52659 | 626 | etaCut = GetEtaSigma(mombin); |
3a385c49 | 627 | } |
3a385c49 | 628 | |
c7d52659 | 629 | if ((phidiff < phiCuthi && phidiff > phiCutlo) && TMath::Abs(etadiff) < etaCut) { |
8c0d179d | 630 | if (track->GetLabel() > fMinMCLabel) { |
631 | ++NMCmatches; | |
632 | trkPMCfrac += track->P(); | |
633 | } | |
875b35fc | 634 | ++Nmatches; |
635 | totalTrkP += track->P(); | |
636 | ||
637 | if (fCreateHisto) { | |
b8b1b9e1 | 638 | if (fHadCorr > 1) { |
875b35fc | 639 | fHistMatchdRvsEP[fCentBin]->Fill(dR, energyclus / mom); |
640 | } | |
1f6c8e39 | 641 | } |
642 | } | |
643 | } | |
8c0d179d | 644 | |
645 | if (totalTrkP > 0) | |
646 | trkPMCfrac /= totalTrkP; | |
3a385c49 | 647 | } |
648 | ||
0b777a09 | 649 | //________________________________________________________________________ |
1f6c8e39 | 650 | Bool_t AliHadCorrTask::Run() |
0b777a09 | 651 | { |
1f6c8e39 | 652 | // Run the hadronic correction |
d6f334b5 | 653 | |
1b4db73a | 654 | if (fCreateHisto) |
655 | DoTrackLoop(); | |
656 | ||
d6f334b5 | 657 | // delete output |
658 | fOutClusters->Delete(); | |
0b777a09 | 659 | |
85835c5b | 660 | // loop over all clusters |
661 | const Int_t Nclus = fCaloClusters->GetEntries(); | |
a89e2005 | 662 | for (Int_t iClus = 0, clusCount=0; iClus < Nclus; ++iClus) { |
4f833c06 | 663 | AliEmcalParticle *emccluster = static_cast<AliEmcalParticle*>(fCaloClusters->At(iClus)); |
1f6c8e39 | 664 | if (!emccluster) |
d6f334b5 | 665 | continue; |
666 | ||
1f6c8e39 | 667 | AliVCluster *cluster = emccluster->GetCluster(); |
668 | if (!cluster) | |
669 | continue; | |
670 | if (!AcceptCluster(cluster)) | |
d6f334b5 | 671 | continue; |
672 | ||
1f6c8e39 | 673 | Double_t energyclus = 0; |
8c0d179d | 674 | if (fCreateHisto) { |
62206c47 | 675 | fHistEbefore->Fill(fCent, cluster->E()); |
8c0d179d | 676 | fHistNclusvsCent->Fill(fCent); |
677 | } | |
32bf39af | 678 | |
d6f334b5 | 679 | // apply correction / subtraction |
0c5f6f08 | 680 | // to subtract only the closest track set fHadCor to a % |
681 | // to subtract all tracks within the cut set fHadCor to %+1 | |
682 | if (fHadCorr > 1) | |
683 | energyclus = ApplyHadCorrAllTracks(emccluster, fHadCorr - 1); | |
f9b8b4b8 | 684 | else if (fHadCorr > 0) |
0c5f6f08 | 685 | energyclus = ApplyHadCorrOneTrack(emccluster, fHadCorr); |
686 | else | |
687 | energyclus = cluster->E(); | |
bcf52e4a | 688 | |
689 | if (energyclus < 0) | |
d6f334b5 | 690 | energyclus = 0; |
d6f334b5 | 691 | |
bcf52e4a | 692 | if (energyclus > 0) { // create corrected cluster |
d6f334b5 | 693 | AliVCluster *oc; |
1b4db73a | 694 | if (fEsdMode) { |
4f833c06 | 695 | AliESDCaloCluster *ec = dynamic_cast<AliESDCaloCluster*>(cluster); |
696 | if (!ec) | |
697 | continue; | |
698 | oc = new ((*fOutClusters)[clusCount]) AliESDCaloCluster(*ec); | |
699 | } else { | |
700 | AliAODCaloCluster *ac = dynamic_cast<AliAODCaloCluster*>(cluster); | |
701 | if (!ac) | |
702 | continue; | |
703 | oc = new ((*fOutClusters)[clusCount]) AliAODCaloCluster(*ac); | |
704 | } | |
d6f334b5 | 705 | oc->SetE(energyclus); |
1f6c8e39 | 706 | |
d6f334b5 | 707 | ++clusCount; |
1f6c8e39 | 708 | |
709 | if (fCreateHisto) | |
710 | fHistEafter->Fill(fCent, energyclus); | |
d6f334b5 | 711 | } |
0b777a09 | 712 | } |
1f6c8e39 | 713 | |
714 | return kTRUE; | |
715 | } | |
bcf52e4a | 716 | |
1f6c8e39 | 717 | //________________________________________________________________________ |
718 | Double_t AliHadCorrTask::ApplyHadCorrOneTrack(AliEmcalParticle *emccluster, Double_t hadCorr) | |
719 | { | |
4f833c06 | 720 | // Apply the hadronic correction with one track only. |
721 | ||
1f6c8e39 | 722 | AliVCluster *cluster = emccluster->GetCluster(); |
181a0a09 | 723 | Double_t energyclus = cluster->E(); |
724 | Int_t iMin = emccluster->GetMatchedObjId(); | |
1f6c8e39 | 725 | if (iMin < 0) |
726 | return energyclus; | |
727 | ||
4f833c06 | 728 | AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iMin)); |
1f6c8e39 | 729 | if (!emctrack) |
730 | return energyclus; | |
731 | ||
181a0a09 | 732 | // check if track also points to cluster |
733 | Int_t cid = emctrack->GetMatchedObjId(); | |
875b35fc | 734 | if (fDoTrackClus && (cid!=emccluster->IdInCollection())) |
c7d52659 | 735 | return energyclus; |
736 | ||
1f6c8e39 | 737 | AliVTrack *track = emctrack->GetTrack(); |
738 | if (!track) | |
739 | return energyclus; | |
740 | ||
741 | Double_t mom = track->P(); | |
742 | if (mom == 0) | |
743 | return energyclus; | |
744 | ||
745 | Double_t dRmin = emccluster->GetMatchedObjDistance(); | |
746 | Double_t dEtaMin = 1e9; | |
747 | Double_t dPhiMin = 1e9; | |
748 | AliPicoTrack::GetEtaPhiDiff(track, cluster, dPhiMin, dEtaMin); | |
749 | ||
b8b1b9e1 | 750 | UInt_t mombin = GetMomBin(mom); |
1f6c8e39 | 751 | Int_t centbinch = fCentBin; |
c7d52659 | 752 | if (track->Charge()<0) |
8c0d179d | 753 | centbinch += fNcentBins; |
1f6c8e39 | 754 | |
181a0a09 | 755 | // plot some histograms if switched on |
1f6c8e39 | 756 | if (fCreateHisto) { |
757 | Int_t etabin = 0; | |
758 | if(track->Eta() > 0) | |
759 | etabin = 1; | |
760 | ||
761 | fHistMatchEtaPhi[centbinch][mombin][etabin]->Fill(dEtaMin, dPhiMin); | |
762 | ||
763 | if (mom > 0) { | |
764 | fHistMatchEvsP[fCentBin]->Fill(energyclus, energyclus / mom); | |
765 | fHistEoPCent->Fill(fCent, energyclus / mom); | |
766 | fHistMatchdRvsEP[fCentBin]->Fill(dRmin, energyclus / mom); | |
767 | } | |
768 | } | |
769 | ||
181a0a09 | 770 | // define eta/phi cuts |
c7d52659 | 771 | Double_t etaCut = 0.0; |
772 | Double_t phiCutlo = 0.0; | |
773 | Double_t phiCuthi = 0.0; | |
1f6c8e39 | 774 | if (fPhiMatch > 0) { |
c7d52659 | 775 | phiCutlo = -fPhiMatch; |
875b35fc | 776 | phiCuthi = +fPhiMatch; |
1f6c8e39 | 777 | } |
778 | else { | |
c7d52659 | 779 | phiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, fCentBin); |
780 | phiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, fCentBin); | |
1f6c8e39 | 781 | } |
782 | if(fEtaMatch > 0) { | |
c7d52659 | 783 | etaCut = fEtaMatch; |
1f6c8e39 | 784 | } |
785 | else { | |
c7d52659 | 786 | etaCut = GetEtaSigma(mombin); |
1f6c8e39 | 787 | } |
788 | ||
181a0a09 | 789 | // apply the correction if the track is in the eta/phi window |
c7d52659 | 790 | if ((dPhiMin < phiCuthi && dPhiMin > phiCutlo) && TMath::Abs(dEtaMin) < etaCut) { |
875b35fc | 791 | energyclus -= hadCorr * mom; |
1f6c8e39 | 792 | } |
793 | ||
794 | return energyclus; | |
0b777a09 | 795 | } |
4711c521 | 796 | |
0b777a09 | 797 | //________________________________________________________________________ |
1f6c8e39 | 798 | Double_t AliHadCorrTask::ApplyHadCorrAllTracks(AliEmcalParticle *emccluster, Double_t hadCorr) |
0b777a09 | 799 | { |
4f833c06 | 800 | // Apply the hadronic correction with all tracks. |
801 | ||
1f6c8e39 | 802 | AliVCluster *cluster = emccluster->GetCluster(); |
803 | ||
804 | Double_t energyclus = cluster->E(); | |
805 | Double_t cNcells = cluster->GetNCells(); | |
806 | ||
807 | Double_t totalTrkP = 0.0; // count total track momentum | |
808 | Int_t Nmatches = 0; // count total number of matches | |
8c0d179d | 809 | |
810 | Double_t trkPMCfrac = 0.0; // count total track momentum | |
811 | Int_t NMCmatches = 0; // count total number of matches | |
1f6c8e39 | 812 | |
181a0a09 | 813 | // do the loop over the matched tracks and get the number of matches and the total momentum |
8c0d179d | 814 | DoMatchedTracksLoop(emccluster, totalTrkP, Nmatches, trkPMCfrac, NMCmatches); |
1f6c8e39 | 815 | |
816 | Double_t Esub = hadCorr * totalTrkP; | |
62206c47 | 817 | |
818 | if (Esub > energyclus) | |
819 | Esub = energyclus; | |
820 | ||
181a0a09 | 821 | // applying Peter's proposed algorithm |
822 | // never subtract the full energy of the cluster | |
1b4db73a | 823 | Double_t clusEexcl = fEexclCell * cNcells; |
824 | if (energyclus < clusEexcl) | |
825 | clusEexcl = energyclus; | |
62206c47 | 826 | if ((energyclus - Esub) < clusEexcl) |
827 | Esub = (energyclus - clusEexcl); | |
828 | ||
1eafa82a | 829 | // embedding |
830 | Double_t EsubMC = 0; | |
831 | Double_t EsubBkg = 0; | |
832 | Double_t EclusMC = 0; | |
833 | Double_t EclusBkg = 0; | |
834 | Double_t EclusCorr = 0; | |
835 | Double_t EclusMCcorr = 0; | |
836 | Double_t EclusBkgcorr = 0; | |
837 | Double_t overSub = 0; | |
838 | if (fIsEmbedded) { | |
839 | EsubMC = hadCorr * totalTrkP * trkPMCfrac; | |
840 | EsubBkg = hadCorr * totalTrkP - EsubMC; | |
841 | EclusMC = energyclus * cluster->GetMCEnergyFraction(); | |
842 | EclusBkg = energyclus - EclusMC; | |
843 | ||
844 | if (energyclus > Esub) | |
845 | EclusCorr = energyclus - Esub; | |
846 | ||
847 | if (EclusMC > EsubMC) | |
848 | EclusMCcorr = EclusMC - EsubMC; | |
849 | ||
850 | if (EclusBkg > EsubBkg) | |
851 | EclusBkgcorr = EclusBkg - EsubBkg; | |
852 | ||
853 | overSub = EclusMCcorr + EclusBkgcorr - EclusCorr; | |
854 | } | |
855 | ||
181a0a09 | 856 | // plot some histograms if switched on |
857 | if (fCreateHisto) { | |
1f6c8e39 | 858 | fHistNMatchCent->Fill(fCent, Nmatches); |
859 | fHistNMatchEnergy[fCentBin]->Fill(energyclus, Nmatches); | |
860 | ||
861 | if(Nmatches > 0) | |
862 | fHistNclusMatchvsCent->Fill(fCent); | |
1b4db73a | 863 | |
864 | if(Nmatches < 3) | |
865 | fHistNCellsEnergy[fCentBin][Nmatches]->Fill(energyclus, cNcells); | |
866 | else | |
1f6c8e39 | 867 | fHistNCellsEnergy[fCentBin][3]->Fill(energyclus, cNcells); |
8c0d179d | 868 | |
869 | if (totalTrkP>0) { | |
870 | Double_t EoP = energyclus / totalTrkP; | |
1f6c8e39 | 871 | fHistEoPCent->Fill(fCent, EoP); |
872 | fHistMatchEvsP[fCentBin]->Fill(energyclus, EoP); | |
8c0d179d | 873 | |
874 | fHistEsubPchRatAll[fCentBin]->Fill(totalTrkP, Esub / totalTrkP); | |
875 | ||
876 | if (Nmatches == 1) { | |
877 | Int_t iMin = emccluster->GetMatchedObjId(); | |
878 | AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iMin)); | |
879 | AliVTrack *track = emctrack->GetTrack(); | |
880 | Int_t centbinchm = fCentBin; | |
881 | if (track->Charge()<0) | |
882 | centbinchm += fNcentBins; | |
883 | ||
884 | fHistEsubPchRat[centbinchm]->Fill(totalTrkP, Esub / totalTrkP); | |
885 | fHistEsubPch[centbinchm]->Fill(totalTrkP, Esub); | |
886 | } | |
1f6c8e39 | 887 | |
1eafa82a | 888 | if (fIsEmbedded) { |
88646d1d | 889 | fHistOversub[fCentBin]->Fill(energyclus, overSub / energyclus); |
890 | ||
891 | if (cluster->GetMCEnergyFraction() > 0.95) | |
892 | fHistOversubMCClusters[fCentBin]->Fill(energyclus, overSub / energyclus); | |
893 | else if (cluster->GetMCEnergyFraction() < 0.05) | |
894 | fHistOversubNonMCClusters[fCentBin]->Fill(energyclus, overSub / energyclus); | |
1eafa82a | 895 | |
896 | if (trkPMCfrac < 0.05) | |
88646d1d | 897 | fHistNonEmbTrackMatchesOversub[fCentBin]->Fill(energyclus, overSub / energyclus); |
1eafa82a | 898 | else if (trkPMCfrac > 0.95) |
88646d1d | 899 | fHistEmbTrackMatchesOversub[fCentBin]->Fill(energyclus, overSub / energyclus); |
56bd3193 | 900 | } |
0c5f6f08 | 901 | } |
1f6c8e39 | 902 | } |
1f6c8e39 | 903 | |
1eafa82a | 904 | if (fIsEmbedded && fDoExact) { |
905 | Esub -= overSub; | |
906 | if (EclusBkgcorr + EclusMCcorr > 0) { | |
907 | Double_t newfrac = EclusMCcorr / (EclusBkgcorr + EclusMCcorr); | |
908 | cluster->SetMCEnergyFraction(newfrac); | |
909 | } | |
910 | } | |
911 | ||
181a0a09 | 912 | // apply the correction |
1f6c8e39 | 913 | energyclus -= Esub; |
914 | ||
915 | return energyclus; | |
0b777a09 | 916 | } |
917 |