Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

mitigate high CPU cost of HBHE hit position calculation in EgammaHcalIsolation #35955

Merged
merged 2 commits into from
Nov 3, 2021
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
Original file line number Diff line number Diff line change
Expand Up @@ -135,7 +135,8 @@ class EgammaHcalIsolation {
}

private:
double goodHitEnergy(const GlobalPoint &pclu,
double goodHitEnergy(float pcluEta,
float pcluPhi,
const HBHERecHit &hit,
int depth,
int ieta,
Expand Down
48 changes: 29 additions & 19 deletions RecoEgamma/EgammaIsolationAlgos/src/EgammaHcalIsolation.cc
Original file line number Diff line number Diff line change
Expand Up @@ -118,52 +118,60 @@ EgammaHcalIsolation::EgammaHcalIsolation(InclusionRule extIncRule,
}
}

double EgammaHcalIsolation::goodHitEnergy(const GlobalPoint &pclu,
double EgammaHcalIsolation::goodHitEnergy(float pcluEta,
float pcluPhi,
const HBHERecHit &hit,
int depth,
int ieta,
int iphi,
int include_or_exclude,
double (*scale)(const double &)) const {
const auto phit = caloGeometry_.getPosition(hit.detid());

if ((extIncRule_ == InclusionRule::withinConeAroundCluster and deltaR2(pclu, phit) > extRadius_) or
(intIncRule_ == InclusionRule::withinConeAroundCluster and deltaR2(pclu, phit) < intRadius_))
return 0.;

const HcalDetId hid(hit.detid());
const int hd = hid.depth(), he = hid.ieta(), hp = hid.iphi();
const int h1 = hd - 1;

if ((hid.subdet() == HcalBarrel and (hd < 1 or hd > int(eThresHB_.size()))) or
(hid.subdet() == HcalEndcap and (hd < 1 or hd > int(eThresHE_.size()))))
edm::LogWarning("EgammaHcalIsolation")
<< " hit in subdet " << hid.subdet() << " has an unaccounted for depth of " << hd << "!!";

if (include_or_exclude == -1 and (he != ieta or hp != iphi))
return 0.;

if (include_or_exclude == 1 and (he == ieta and hp == iphi))
return 0.;

if ((hid.subdet() == HcalBarrel and (hd < 1 or hd > int(eThresHB_.size()))) or
(hid.subdet() == HcalEndcap and (hd < 1 or hd > int(eThresHE_.size()))))
edm::LogWarning("EgammaHcalIsolation")
<< " hit in subdet " << hid.subdet() << " has an unaccounted for depth of " << hd << "!!";

const bool right_depth = (depth == 0 or hd == depth);
if (!right_depth)
return 0.;

const bool goodHBe = hid.subdet() == HcalBarrel and hit.energy() > eThresHB_[h1];
const bool goodHEe = hid.subdet() == HcalEndcap and hit.energy() > eThresHE_[h1];
if (!(goodHBe or goodHEe))
return 0.;

const auto phit = caloGeometry_.getPosition(hit.detid());
Copy link
Contributor

@VinInn VinInn Nov 12, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

suggested change:

auto const & phit = caloGeometry_.getGeometry(hit.detid())->repPos();

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

suggested change:

Thanks. Testing.
PR soon.

const float phitEta = phit.eta();
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

use cashed eta in caloGeometry_

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

please clarify the member method
Geometry/CaloGeometry/interface/CaloGeometry.h has only xyz position

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

RhoEtaPhi const& repPos() const { return m_rep; }
virtual float rhoPos() const { return m_rep.rho(); }
virtual float etaPos() const { return m_rep.eta(); }
virtual float phiPos() const { return m_rep.phi(); }

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

eventually through CaloGeometry::getGeometry if not in the same SubDetector

Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

one could extend the CaloGeometry interface to support repPos (copy paste code of getPosition)
In my opinion would be more efficient to check first if the CellGeometry is the same and avoid the overhead of the code (not inlined) to search for the subDetector for each single hit....

Copy link
Contributor Author

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

ah, so the suggestion is for an extension to the geometry interface; I thought at first that this is something that exists already.

IIRC, there is one rechit per cell in ECAL/HCAL

Copy link
Contributor

@VinInn VinInn Nov 12, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

it exits already in CaloCellGeometry that is what CaloGeometry access. please have a look to the implementation in
https://cmssdt.cern.ch/dxr/CMSSW/source/Geometry/CaloGeometry/src/CaloGeometry.cc#50
here locally one can invoke
caloGeometry_.getGeometry(hit.detid())->repPos();

Copy link
Contributor

@VinInn VinInn Nov 12, 2021

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

here https://github.com/slava77/cmssw/blob/f33bdaa8be9c4aa4a3f336dbf23501b059a8988b/RecoEgamma/EgammaIsolationAlgos/src/EgammaHcalIsolation.cc#L188
thre is a loop on the hits in HBorHE. one can verify if a hit is in the same SubDet of the previous one and avoid all the search logic in CaloGeometry and acess directly CaloCellGeometry


if (extIncRule_ == InclusionRule::withinConeAroundCluster or intIncRule_ == InclusionRule::withinConeAroundCluster) {
auto const dR2 = deltaR2(pcluEta, pcluPhi, phitEta, phit.phi());
Copy link
Contributor

Choose a reason for hiding this comment

The reason will be displayed to describe this comment to others. Learn more.

use cashed phi in caloGeometry

if ((extIncRule_ == InclusionRule::withinConeAroundCluster and dR2 > extRadius_) or
(intIncRule_ == InclusionRule::withinConeAroundCluster and dR2 < intRadius_))
return 0.;
}

DetId did = hcalTopology_.idFront(hid);
const uint32_t flag = hit.flags();
const uint32_t dbflag = hcalChStatus_.getValues(did)->getValue();
int severity = hcalSevLvlComputer_.getSeverityLevel(did, flag, dbflag);
bool recovered = hcalSevLvlComputer_.recoveredRecHit(did, flag);

const double het = hit.energy() * scaleToEt(phit.eta());
const bool goodHB = hid.subdet() == HcalBarrel and (severity <= maxSeverityHB_ or recovered) and
hit.energy() > eThresHB_[h1] and het > etThresHB_[h1];
const bool goodHE = hid.subdet() == HcalEndcap and (severity <= maxSeverityHE_ or recovered) and
hit.energy() > eThresHE_[h1] and het > etThresHE_[h1];
const double het = hit.energy() * scaleToEt(phitEta);
const bool goodHB = goodHBe and (severity <= maxSeverityHB_ or recovered) and het > etThresHB_[h1];
const bool goodHE = goodHEe and (severity <= maxSeverityHE_ or recovered) and het > etThresHE_[h1];

if (goodHB or goodHE)
return hit.energy() * scale(phit.eta());
return hit.energy() * scale(phitEta);

return 0.;
}
Expand All @@ -175,8 +183,10 @@ double EgammaHcalIsolation::getHcalSum(const GlobalPoint &pclu,
int include_or_exclude,
double (*scale)(const double &)) const {
double sum = 0.;
const float pcluEta = pclu.eta();
const float pcluPhi = pclu.phi();
for (const auto &hit : mhbhe_)
sum += goodHitEnergy(pclu, hit, depth, ieta, iphi, include_or_exclude, scale);
sum += goodHitEnergy(pcluEta, pcluPhi, hit, depth, ieta, iphi, include_or_exclude, scale);

return sum;
}