3 // Hadronic correction task.
5 // Author: R.Reed, C.Loizides, S.Aiola
8 #include <TClonesArray.h>
13 #include "AliAnalysisManager.h"
14 #include "AliAODEvent.h"
15 #include "AliAODCaloCluster.h"
16 #include "AliESDCaloCluster.h"
17 #include "AliVTrack.h"
18 #include "AliPicoTrack.h"
19 #include "AliVEventHandler.h"
20 #include "AliEmcalParticle.h"
21 #include "AliEMCALGeometry.h"
23 #include "AliHadCorrTask.h"
25 ClassImp(AliHadCorrTask)
27 //________________________________________________________________________
28 AliHadCorrTask::AliHadCorrTask() :
29 AliAnalysisTaskEmcal("AliHadCorrTask", kFALSE),
40 fHistNclusMatchvsCent(0),
45 fHistNClusMatchCent(0)
47 // Default constructor.
49 for(Int_t i=0; i<8; i++) {
51 fHistEsubPchRat[i] = 0;
52 fHistEsubPchRatAll[i] = 0;
55 fHistMatchEvsP[i] = 0;
56 fHistMatchdRvsEP[i] = 0;
57 fHistNMatchEnergy[i] = 0;
59 fHistEmbTrackMatchesOversub[i] = 0;
60 fHistNonEmbTrackMatchesOversub[i] = 0;
61 fHistOversubMCClusters[i] = 0;
62 fHistOversubNonMCClusters[i] = 0;
65 for(Int_t j=0; j<4; j++)
66 fHistNCellsEnergy[i][j] = 0;
69 for(Int_t j=0; j<9; j++) {
70 for(Int_t k=0; k<2; k++)
71 fHistMatchEtaPhi[i][j][k] = 0;
76 //________________________________________________________________________
77 AliHadCorrTask::AliHadCorrTask(const char *name, Bool_t histo) :
78 AliAnalysisTaskEmcal(name, histo),
79 fOutCaloName("CaloClustersCorr"),
89 fHistNclusMatchvsCent(0),
94 fHistNClusMatchCent(0)
96 // Standard constructor.
98 for(Int_t i=0; i<8; i++) {
100 fHistEsubPchRat[i] = 0;
101 fHistEsubPchRatAll[i] = 0;
104 fHistMatchEvsP[i] = 0;
105 fHistMatchdRvsEP[i] = 0;
106 fHistNMatchEnergy[i] = 0;
108 fHistEmbTrackMatchesOversub[i] = 0;
109 fHistNonEmbTrackMatchesOversub[i] = 0;
110 fHistOversubMCClusters[i] = 0;
111 fHistOversubNonMCClusters[i] = 0;
114 for(Int_t j=0; j<4; j++)
115 fHistNCellsEnergy[i][j] = 0;
118 for(Int_t j=0; j<9; j++) {
119 for(Int_t k=0; k<2; k++)
120 fHistMatchEtaPhi[i][j][k] = 0;
124 SetMakeGeneralHistograms(histo);
126 fBranchNames="ESD:AliESDRun.,AliESDHeader.,PrimaryVertex.";
129 //________________________________________________________________________
130 AliHadCorrTask::~AliHadCorrTask()
135 //________________________________________________________________________
136 UInt_t AliHadCorrTask::GetMomBin(Double_t p) const
143 else if (p>=0.5 && p<1.0)
145 else if (p>=1.0 && p<1.5)
147 else if (p>=1.5 && p<2.)
149 else if (p>=2. && p<3.)
151 else if (p>=3. && p<4.)
153 else if (p>=4. && p<5.)
155 else if (p>=5. && p<8.)
163 //________________________________________________________________________
164 Double_t AliHadCorrTask::GetEtaSigma(Int_t pbin) const
168 Double_t EtaSigma[9]={0.0097,0.0075,0.0059,0.0055,0.0053,0.005,0.005,0.0045,0.0042};
169 return 2.0*EtaSigma[pbin];
172 //________________________________________________________________________
173 Double_t AliHadCorrTask::GetPhiMean(Int_t pbin, Int_t centbin) const
178 Double_t PhiMean[9]={0.036,
187 return PhiMean[pbin];
188 } else if (centbin==1) {
189 Double_t PhiMean[9]={0.036,
198 return PhiMean[pbin];
199 } else if (centbin==2) {
200 Double_t PhiMean[9]={0.036,
209 return PhiMean[pbin];
210 } else if (centbin==3) {
211 Double_t PhiMean[9]={0.036,
220 return PhiMean[pbin];
221 } else if (centbin==4) {
222 Double_t PhiMean[9]={0.036,
231 return PhiMean[pbin]*(-1.);
232 } else if (centbin==5) {
233 Double_t PhiMean[9]={0.036,
242 return PhiMean[pbin]*(-1.);
243 } else if (centbin==6) {
244 Double_t PhiMean[9]={0.036,
253 return PhiMean[pbin]*(-1.);
254 } else if (centbin==7) {
255 Double_t PhiMean[9]={0.036,
264 return PhiMean[pbin]*(-1.);
270 //________________________________________________________________________
271 Double_t AliHadCorrTask::GetPhiSigma(Int_t pbin, Int_t centbin) const
276 Double_t PhiSigma[9]={0.0221,
285 return 2.*PhiSigma[pbin];
286 } else if (centbin==1) {
287 Double_t PhiSigma[9]={0.0217,
296 return 2.*PhiSigma[pbin];
297 } else if (centbin==2) {
298 Double_t PhiSigma[9]={0.0211,
307 return 2.*PhiSigma[pbin];
308 } else if (centbin==3) {
309 Double_t PhiSigma[9]={0.0215,
318 return 2.*PhiSigma[pbin];
319 } else if (centbin==4) {
320 Double_t PhiSigma[9]={0.0199,
329 return 2.*PhiSigma[pbin];
330 } else if (centbin==5) {
331 Double_t PhiSigma[9]={0.0200,
340 return 2.*PhiSigma[pbin];
341 } else if (centbin==6) {
342 Double_t PhiSigma[9]={0.0202,
351 return 2.*PhiSigma[pbin];
352 } else if (centbin==7) {
353 Double_t PhiSigma[9]={0.0205,
362 return 2.*PhiSigma[pbin];
368 //________________________________________________________________________
369 void AliHadCorrTask::UserCreateOutputObjects()
371 // Create my user objects.
376 AliAnalysisTaskEmcal::UserCreateOutputObjects();
380 const Int_t nCentChBins = fNcentBins * 2;
382 for(Int_t icent=0; icent<nCentChBins; ++icent) {
383 for(Int_t ipt=0; ipt<9; ++ipt) {
384 for(Int_t ieta=0; ieta<2; ++ieta) {
385 name = Form("fHistMatchEtaPhi_%i_%i_%i",icent,ipt,ieta);
386 fHistMatchEtaPhi[icent][ipt][ieta] = new TH2F(name, name, fNbins, -0.1, 0.1, fNbins, -0.1, 0.1);
387 fOutput->Add(fHistMatchEtaPhi[icent][ipt][ieta]);
391 name = Form("fHistEsubPch_%i",icent);
392 fHistEsubPch[icent]=new TH1F(name, name, fNbins, fMinBinPt, fMaxBinPt);
393 fOutput->Add(fHistEsubPch[icent]);
395 name = Form("fHistEsubPchRat_%i",icent);
396 fHistEsubPchRat[icent]=new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.);
397 fOutput->Add(fHistEsubPchRat[icent]);
399 name = Form("fHistEsubPchRatAll_%i",icent);
400 fHistEsubPchRatAll[icent]=new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.);
401 fOutput->Add(fHistEsubPchRatAll[icent]);
403 if (icent<fNcentBins) {
404 for(Int_t itrk=0; itrk<4; ++itrk) {
405 name = Form("fHistNCellsEnergy_%i_%i",icent,itrk);
406 fHistNCellsEnergy[icent][itrk] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, 101, -0.5, 100.5);
407 fOutput->Add(fHistNCellsEnergy[icent][itrk]);
410 name = Form("fHistMatchEvsP_%i",icent);
411 fHistMatchEvsP[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins*2, 0., 10.);
412 fOutput->Add(fHistMatchEvsP[icent]);
414 name = Form("fHistMatchdRvsEP_%i",icent);
415 fHistMatchdRvsEP[icent] = new TH2F(name, name, fNbins, 0., 0.2, fNbins*2, 0., 10.);
416 fOutput->Add(fHistMatchdRvsEP[icent]);
418 name = Form("fHistNMatchEnergy_%i",icent);
419 fHistNMatchEnergy[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, 101, -0.5, 100.5);
420 fOutput->Add(fHistNMatchEnergy[icent]);
424 fHistNclusvsCent = new TH1F("Nclusvscent", "NclusVsCent", 100, 0, 100);
425 fHistNclusMatchvsCent = new TH1F("NclusMatchvscent", "NclusMatchVsCent", 100, 0, 100);
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);
432 fOutput->Add(fHistNclusMatchvsCent);
433 fOutput->Add(fHistNclusvsCent);
434 fOutput->Add(fHistEbefore);
435 fOutput->Add(fHistEafter);
436 fOutput->Add(fHistEoPCent);
437 fOutput->Add(fHistNMatchCent);
438 fOutput->Add(fHistNClusMatchCent);
441 for(Int_t icent=0; icent<fNcentBins; ++icent) {
442 name = Form("fHistEmbTrackMatchesOversub_%d",icent);
443 fHistEmbTrackMatchesOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
444 fHistEmbTrackMatchesOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
445 fHistEmbTrackMatchesOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
446 fOutput->Add(fHistEmbTrackMatchesOversub[icent]);
448 name = Form("fHistNonEmbTrackMatchesOversub_%d",icent);
449 fHistNonEmbTrackMatchesOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
450 fHistNonEmbTrackMatchesOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
451 fHistNonEmbTrackMatchesOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
452 fOutput->Add(fHistNonEmbTrackMatchesOversub[icent]);
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]);
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]);
466 name = Form("fHistOversub_%d",icent);
467 fHistOversub[icent] = new TH2F(name, name, fNbins, fMinBinPt, fMaxBinPt, fNbins, 0, 1.2);
468 fHistOversub[icent]->GetXaxis()->SetTitle("E_{clus}^{raw} (GeV)");
469 fHistOversub[icent]->GetYaxis()->SetTitle("E_{oversub} / E_{clus}^{raw}");
470 fOutput->Add(fHistOversub[icent]);
474 PostData(1, fOutput);
477 //________________________________________________________________________
478 void AliHadCorrTask::ExecOnce()
480 // Do base class initializations and if it fails -> bail out
482 AliAnalysisTaskEmcal::ExecOnce();
487 if (dynamic_cast<AliAODEvent*>(InputEvent()))
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.");
500 fOutClusters = new TClonesArray("AliESDCaloCluster");
502 fOutClusters = new TClonesArray("AliAODCaloCluster");
504 fOutClusters->SetName(fOutCaloName);
506 // post output in event if not yet present
507 if (!(InputEvent()->FindListObject(fOutCaloName))) {
508 InputEvent()->AddObject(fOutClusters);
511 fInitialized = kFALSE;
512 AliFatal(Form("%s: Container with same name %s already present. Aborting", GetName(), fOutCaloName.Data()));
517 //________________________________________________________________________
518 void AliHadCorrTask::DoTrackLoop()
520 const Int_t Ntracks = fTracks->GetEntries();
521 for (Int_t iTrk = 0; iTrk < Ntracks; ++iTrk) {
522 Int_t NmatchClus = 0;
524 AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iTrk));
527 if (!emctrack->IsEMCAL())
530 AliVTrack *track = emctrack->GetTrack();
533 if (!AcceptTrack(track))
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)
542 Int_t Nclus = emctrack->GetNumberOfMatchedObj();
544 for (Int_t iClus = 0; iClus < Nclus; ++iClus) {
545 AliEmcalParticle *emccluster = static_cast<AliEmcalParticle*>(fCaloClusters->At(iClus));
549 AliVCluster *cluster = emccluster->GetCluster();
552 if (!AcceptCluster(cluster))
555 Double_t etadiff = 999;
556 Double_t phidiff = 999;
557 AliPicoTrack::GetEtaPhiDiff(track, cluster, phidiff, etadiff);
559 if (TMath::Abs(phidiff) < fPhiMatch && TMath::Abs(etadiff) < fEtaMatch)
562 fHistNClusMatchCent->Fill(fCent, NmatchClus);
566 //________________________________________________________________________
567 void AliHadCorrTask::DoMatchedTracksLoop(AliEmcalParticle *emccluster, Double_t &totalTrkP, Int_t &Nmatches, Double_t &trkPMCfrac, Int_t &NMCmatches)
569 // Do the loop over matched tracks for cluster emccluster.
571 AliVCluster *cluster = emccluster->GetCluster();
572 Int_t iClus = emccluster->IdInCollection();
573 Double_t energyclus = cluster->E();
575 // loop over matched tracks
576 const Int_t Ntrks = emccluster->GetNumberOfMatchedObj();
577 for (Int_t i = 0; i < Ntrks; ++i) {
578 Int_t iTrack = emccluster->GetMatchedObjId(i);
579 Double_t dR = emccluster->GetMatchedObjDistance(i);
581 AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iTrack));
585 // check if track also points to cluster
586 if (fDoTrackClus && (emctrack->GetMatchedObjId(0)) != iClus)
589 AliVTrack *track = emctrack->GetTrack();
592 if (!AcceptTrack(track))
595 Double_t etadiff = 999;
596 Double_t phidiff = 999;
597 AliPicoTrack::GetEtaPhiDiff(track, cluster, phidiff, etadiff);
599 Double_t mom = track->P();
600 UInt_t mombin = GetMomBin(mom);
601 Int_t centbinch = fCentBin;
602 if (track->Charge()<0)
603 centbinch += fNcentBins;
607 if(track->Eta() > 0) etabin=1;
608 fHistMatchEtaPhi[centbinch][mombin][etabin]->Fill(etadiff, phidiff);
611 Double_t etaCut = 0.0;
612 Double_t phiCutlo = 0.0;
613 Double_t phiCuthi = 0.0;
616 phiCutlo = -fPhiMatch;
617 phiCuthi = +fPhiMatch;
619 phiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, fCentBin);
620 phiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, fCentBin);
626 etaCut = GetEtaSigma(mombin);
629 if ((phidiff < phiCuthi && phidiff > phiCutlo) && TMath::Abs(etadiff) < etaCut) {
630 if (track->GetLabel() > fMinMCLabel) {
632 trkPMCfrac += track->P();
635 totalTrkP += track->P();
639 fHistMatchdRvsEP[fCentBin]->Fill(dR, energyclus / mom);
646 trkPMCfrac /= totalTrkP;
649 //________________________________________________________________________
650 Bool_t AliHadCorrTask::Run()
652 // Run the hadronic correction
658 fOutClusters->Delete();
660 // loop over all clusters
661 const Int_t Nclus = fCaloClusters->GetEntries();
662 for (Int_t iClus = 0, clusCount=0; iClus < Nclus; ++iClus) {
663 AliEmcalParticle *emccluster = static_cast<AliEmcalParticle*>(fCaloClusters->At(iClus));
667 AliVCluster *cluster = emccluster->GetCluster();
670 if (!AcceptCluster(cluster))
673 Double_t energyclus = 0;
675 fHistEbefore->Fill(fCent, cluster->E());
676 fHistNclusvsCent->Fill(fCent);
679 // apply correction / subtraction
680 // to subtract only the closest track set fHadCor to a %
681 // to subtract all tracks within the cut set fHadCor to %+1
683 energyclus = ApplyHadCorrAllTracks(emccluster, fHadCorr - 1);
684 else if (fHadCorr > 0)
685 energyclus = ApplyHadCorrOneTrack(emccluster, fHadCorr);
687 energyclus = cluster->E();
692 if (energyclus > 0) { // create corrected cluster
695 AliESDCaloCluster *ec = dynamic_cast<AliESDCaloCluster*>(cluster);
698 oc = new ((*fOutClusters)[clusCount]) AliESDCaloCluster(*ec);
700 AliAODCaloCluster *ac = dynamic_cast<AliAODCaloCluster*>(cluster);
703 oc = new ((*fOutClusters)[clusCount]) AliAODCaloCluster(*ac);
705 oc->SetE(energyclus);
710 fHistEafter->Fill(fCent, energyclus);
717 //________________________________________________________________________
718 Double_t AliHadCorrTask::ApplyHadCorrOneTrack(AliEmcalParticle *emccluster, Double_t hadCorr)
720 // Apply the hadronic correction with one track only.
722 AliVCluster *cluster = emccluster->GetCluster();
723 Double_t energyclus = cluster->E();
724 Int_t iMin = emccluster->GetMatchedObjId();
728 AliEmcalParticle *emctrack = static_cast<AliEmcalParticle*>(fTracks->At(iMin));
732 // check if track also points to cluster
733 Int_t cid = emctrack->GetMatchedObjId();
734 if (fDoTrackClus && (cid!=emccluster->IdInCollection()))
737 AliVTrack *track = emctrack->GetTrack();
741 Double_t mom = track->P();
745 Double_t dRmin = emccluster->GetMatchedObjDistance();
746 Double_t dEtaMin = 1e9;
747 Double_t dPhiMin = 1e9;
748 AliPicoTrack::GetEtaPhiDiff(track, cluster, dPhiMin, dEtaMin);
750 UInt_t mombin = GetMomBin(mom);
751 Int_t centbinch = fCentBin;
752 if (track->Charge()<0)
753 centbinch += fNcentBins;
755 // plot some histograms if switched on
761 fHistMatchEtaPhi[centbinch][mombin][etabin]->Fill(dEtaMin, dPhiMin);
764 fHistMatchEvsP[fCentBin]->Fill(energyclus, energyclus / mom);
765 fHistEoPCent->Fill(fCent, energyclus / mom);
766 fHistMatchdRvsEP[fCentBin]->Fill(dRmin, energyclus / mom);
770 // define eta/phi cuts
771 Double_t etaCut = 0.0;
772 Double_t phiCutlo = 0.0;
773 Double_t phiCuthi = 0.0;
775 phiCutlo = -fPhiMatch;
776 phiCuthi = +fPhiMatch;
779 phiCutlo = GetPhiMean(mombin, centbinch) - GetPhiSigma(mombin, fCentBin);
780 phiCuthi = GetPhiMean(mombin, centbinch) + GetPhiSigma(mombin, fCentBin);
786 etaCut = GetEtaSigma(mombin);
789 // apply the correction if the track is in the eta/phi window
790 if ((dPhiMin < phiCuthi && dPhiMin > phiCutlo) && TMath::Abs(dEtaMin) < etaCut) {
791 energyclus -= hadCorr * mom;
797 //________________________________________________________________________
798 Double_t AliHadCorrTask::ApplyHadCorrAllTracks(AliEmcalParticle *emccluster, Double_t hadCorr)
800 // Apply the hadronic correction with all tracks.
802 AliVCluster *cluster = emccluster->GetCluster();
804 Double_t energyclus = cluster->E();
805 Double_t cNcells = cluster->GetNCells();
807 Double_t totalTrkP = 0.0; // count total track momentum
808 Int_t Nmatches = 0; // count total number of matches
810 Double_t trkPMCfrac = 0.0; // count total track momentum
811 Int_t NMCmatches = 0; // count total number of matches
813 // do the loop over the matched tracks and get the number of matches and the total momentum
814 DoMatchedTracksLoop(emccluster, totalTrkP, Nmatches, trkPMCfrac, NMCmatches);
816 Double_t Esub = hadCorr * totalTrkP;
818 if (Esub > energyclus)
821 // applying Peter's proposed algorithm
822 // never subtract the full energy of the cluster
823 Double_t clusEexcl = fEexclCell * cNcells;
824 if (energyclus < clusEexcl)
825 clusEexcl = energyclus;
826 if ((energyclus - Esub) < clusEexcl)
827 Esub = (energyclus - clusEexcl);
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;
839 EsubMC = hadCorr * totalTrkP * trkPMCfrac;
840 EsubBkg = hadCorr * totalTrkP - EsubMC;
841 EclusMC = energyclus * cluster->GetMCEnergyFraction();
842 EclusBkg = energyclus - EclusMC;
844 if (energyclus > Esub)
845 EclusCorr = energyclus - Esub;
847 if (EclusMC > EsubMC)
848 EclusMCcorr = EclusMC - EsubMC;
850 if (EclusBkg > EsubBkg)
851 EclusBkgcorr = EclusBkg - EsubBkg;
853 overSub = EclusMCcorr + EclusBkgcorr - EclusCorr;
856 // plot some histograms if switched on
858 fHistNMatchCent->Fill(fCent, Nmatches);
859 fHistNMatchEnergy[fCentBin]->Fill(energyclus, Nmatches);
862 fHistNclusMatchvsCent->Fill(fCent);
865 fHistNCellsEnergy[fCentBin][Nmatches]->Fill(energyclus, cNcells);
867 fHistNCellsEnergy[fCentBin][3]->Fill(energyclus, cNcells);
870 Double_t EoP = energyclus / totalTrkP;
871 fHistEoPCent->Fill(fCent, EoP);
872 fHistMatchEvsP[fCentBin]->Fill(energyclus, EoP);
874 fHistEsubPchRatAll[fCentBin]->Fill(totalTrkP, Esub / totalTrkP);
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;
884 fHistEsubPchRat[centbinchm]->Fill(totalTrkP, Esub / totalTrkP);
885 fHistEsubPch[centbinchm]->Fill(totalTrkP, Esub);
889 fHistOversub[fCentBin]->Fill(energyclus, overSub / energyclus);
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);
896 if (trkPMCfrac < 0.05)
897 fHistNonEmbTrackMatchesOversub[fCentBin]->Fill(energyclus, overSub / energyclus);
898 else if (trkPMCfrac > 0.95)
899 fHistEmbTrackMatchesOversub[fCentBin]->Fill(energyclus, overSub / energyclus);
904 if (fIsEmbedded && fDoExact) {
906 if (EclusBkgcorr + EclusMCcorr > 0) {
907 Double_t newfrac = EclusMCcorr / (EclusBkgcorr + EclusMCcorr);
908 cluster->SetMCEnergyFraction(newfrac);
912 // apply the correction