]> git.uio.no Git - u/mrichter/AliRoot.git/blame - PWGGA/EMCALJetTasks/AliHadCorrTask.cxx
git://git.uio.no / u / mrichter / AliRoot.git / blame
? search:
summary | shortlog | log | commit | commitdiff | tree
blob | blame (incremental) | history | HEAD
histo switch
[u/mrichter/AliRoot.git] / PWGGA / EMCALJetTasks / AliHadCorrTask.cxx
CommitLineData
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
25ClassImp(AliHadCorrTask)
26
27//________________________________________________________________________
d6f334b5 28AliHadCorrTask::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//________________________________________________________________________
69AliHadCorrTask::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//________________________________________________________________________
112AliHadCorrTask::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//________________________________________________________________________
158AliHadCorrTask::~AliHadCorrTask()
159{
160 // Destructor
161}
162
d6f334b5 163//________________________________________________________________________
164Int_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//________________________________________________________________________
182Int_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//________________________________________________________________________
210Double_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//________________________________________________________________________
217Double_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//________________________________________________________________________
314Double_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,
3370.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//________________________________________________________________________
412void 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//________________________________________________________________________
495void 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) {
532 AliError(Form("Centrality negative: %f", cent));
533 return;
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//________________________________________________________________________
802void AliHadCorrTask::Terminate(Option_t *)
803{
d6f334b5 804 // Nothing to be done.
0b777a09 805}
806