Skip to content

Commit

Permalink
Merge pull request #44693 from jsamudio/fixPFHCALoff
Browse files Browse the repository at this point in the history
Fix crash in PF Clustering due to empty HCAL RecHits collection
  • Loading branch information
cmsbuild authored Apr 11, 2024
2 parents 1471374 + 0398bdc commit fc8422f
Show file tree
Hide file tree
Showing 4 changed files with 67 additions and 48 deletions.
Original file line number Diff line number Diff line change
Expand Up @@ -197,41 +197,45 @@ void LegacyPFClusterProducer::produce(edm::Event& event, const edm::EventSetup&
auto const& pfClusterSoA = event.get(pfClusterSoAToken_).const_view();
auto const& pfRecHitFractionSoA = event.get(pfRecHitFractionSoAToken_).const_view();

int nRH = pfRecHits.view().size();
int nRH = 0;
if (pfRecHits->metadata().size() != 0)
nRH = pfRecHits.view().size();
reco::PFClusterCollection out;
out.reserve(nRH);

auto const rechitsHandle = event.getHandle(recHitsLabel_);

// Build PFClusters in legacy format
std::vector<int> nTopoSeeds(nRH, 0);
if (nRH != 0) {
// Build PFClusters in legacy format
std::vector<int> nTopoSeeds(nRH, 0);

for (int i = 0; i < pfClusterSoA.nSeeds(); i++) {
nTopoSeeds[pfClusterSoA[i].topoId()]++;
}
for (int i = 0; i < pfClusterSoA.nSeeds(); i++) {
nTopoSeeds[pfClusterSoA[i].topoId()]++;
}

// Looping over SoA PFClusters to produce legacy PFCluster collection
for (int i = 0; i < pfClusterSoA.nSeeds(); i++) {
unsigned int n = pfClusterSoA[i].seedRHIdx();
reco::PFCluster temp;
temp.setSeed((*rechitsHandle)[n].detId()); // Pulling the detId of this PFRecHit from the legacy format input
int offset = pfClusterSoA[i].rhfracOffset();
for (int k = offset; k < (offset + pfClusterSoA[i].rhfracSize()) && k >= 0;
k++) { // Looping over PFRecHits in the same topo cluster
if (pfRecHitFractionSoA[k].pfrhIdx() < nRH && pfRecHitFractionSoA[k].pfrhIdx() > -1 &&
pfRecHitFractionSoA[k].frac() > 0.0) {
const reco::PFRecHitRef& refhit = reco::PFRecHitRef(rechitsHandle, pfRecHitFractionSoA[k].pfrhIdx());
temp.addRecHitFraction(reco::PFRecHitFraction(refhit, pfRecHitFractionSoA[k].frac()));
// Looping over SoA PFClusters to produce legacy PFCluster collection
for (int i = 0; i < pfClusterSoA.nSeeds(); i++) {
unsigned int n = pfClusterSoA[i].seedRHIdx();
reco::PFCluster temp;
temp.setSeed((*rechitsHandle)[n].detId()); // Pulling the detId of this PFRecHit from the legacy format input
int offset = pfClusterSoA[i].rhfracOffset();
for (int k = offset; k < (offset + pfClusterSoA[i].rhfracSize()) && k >= 0;
k++) { // Looping over PFRecHits in the same topo cluster
if (pfRecHitFractionSoA[k].pfrhIdx() < nRH && pfRecHitFractionSoA[k].pfrhIdx() > -1 &&
pfRecHitFractionSoA[k].frac() > 0.0) {
const reco::PFRecHitRef& refhit = reco::PFRecHitRef(rechitsHandle, pfRecHitFractionSoA[k].pfrhIdx());
temp.addRecHitFraction(reco::PFRecHitFraction(refhit, pfRecHitFractionSoA[k].frac()));
}
}
}

// Now PFRecHitFraction of this PFCluster is set. Now compute calculateAndSetPosition (energy, position etc)
if (nTopoSeeds[pfClusterSoA[i].topoId()] == 1 && allCellsPositionCalc_) {
allCellsPositionCalc_->calculateAndSetPosition(temp, paramPF);
} else {
positionCalc_->calculateAndSetPosition(temp, paramPF);
// Now PFRecHitFraction of this PFCluster is set. Now compute calculateAndSetPosition (energy, position etc)
if (nTopoSeeds[pfClusterSoA[i].topoId()] == 1 && allCellsPositionCalc_) {
allCellsPositionCalc_->calculateAndSetPosition(temp, paramPF);
} else {
positionCalc_->calculateAndSetPosition(temp, paramPF);
}
out.emplace_back(std::move(temp));
}
out.emplace_back(std::move(temp));
}

event.emplace(legacyPfClustersToken_, std::move(out));
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -29,16 +29,26 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE {
const reco::PFClusterParamsDeviceCollection& params = setup.getData(pfClusParamsToken);
const reco::PFRecHitHCALTopologyDeviceCollection& topology = setup.getData(topologyToken_);
const reco::PFRecHitHostCollection& pfRecHits = event.get(inputPFRecHitSoA_Token_);
const int nRH = pfRecHits->size();
int nRH = 0;
if (pfRecHits->metadata().size() != 0)
nRH = pfRecHits->size();

reco::PFClusteringVarsDeviceCollection pfClusteringVars{nRH, event.queue()};
reco::PFClusteringEdgeVarsDeviceCollection pfClusteringEdgeVars{(nRH * 8), event.queue()};
reco::PFClusterDeviceCollection pfClusters{nRH, event.queue()};
reco::PFRecHitFractionDeviceCollection pfrhFractions{nRH * pfRecHitFractionAllocation_, event.queue()};

PFClusterProducerKernel kernel(event.queue(), pfRecHits);
kernel.execute(
event.queue(), params, topology, pfClusteringVars, pfClusteringEdgeVars, pfRecHits, pfClusters, pfrhFractions);
if (nRH != 0) {
PFClusterProducerKernel kernel(event.queue(), pfRecHits);
kernel.execute(event.queue(),
params,
topology,
pfClusteringVars,
pfClusteringEdgeVars,
pfRecHits,
pfClusters,
pfrhFractions);
}

if (synchronise_)
alpaka::wait(event.queue());
Expand Down
32 changes: 17 additions & 15 deletions RecoParticleFlow/PFRecHitProducer/plugins/LegacyPFRecHitProducer.cc
Original file line number Diff line number Diff line change
Expand Up @@ -54,23 +54,25 @@ void LegacyPFRecHitProducer::produce(edm::Event& event, const edm::EventSetup& s
const reco::PFRecHitHostCollection::ConstView& alpakaPfRecHits = pfRecHitsAlpakaSoA.const_view();

reco::PFRecHitCollection out;
out.reserve(alpakaPfRecHits.size());

for (size_t i = 0; i < alpakaPfRecHits.size(); i++) {
reco::PFRecHit& pfrh =
out.emplace_back(caloGeo_.at(alpakaPfRecHits[i].layer())->getGeometry(alpakaPfRecHits[i].detId()),
alpakaPfRecHits[i].detId(),
alpakaPfRecHits[i].layer(),
alpakaPfRecHits[i].energy());
pfrh.setTime(alpakaPfRecHits[i].time());
pfrh.setDepth(alpakaPfRecHits[i].depth());
if (alpakaPfRecHits.metadata().size() != 0) {
out.reserve(alpakaPfRecHits.size());
for (size_t i = 0; i < alpakaPfRecHits.size(); i++) {
reco::PFRecHit& pfrh =
out.emplace_back(caloGeo_.at(alpakaPfRecHits[i].layer())->getGeometry(alpakaPfRecHits[i].detId()),
alpakaPfRecHits[i].detId(),
alpakaPfRecHits[i].layer(),
alpakaPfRecHits[i].energy());
pfrh.setTime(alpakaPfRecHits[i].time());
pfrh.setDepth(alpakaPfRecHits[i].depth());

// order in Alpaka: N, S, E, W,NE,SW,SE,NW
const short eta[8] = {0, 0, 1, -1, 1, -1, 1, -1};
const short phi[8] = {1, -1, 0, 0, 1, -1, -1, 1};
for (size_t k = 0; k < 8; k++)
if (alpakaPfRecHits[i].neighbours()(k) != -1)
pfrh.addNeighbour(eta[k], phi[k], 0, alpakaPfRecHits[i].neighbours()(k));
// order in Alpaka: N, S, E, W,NE,SW,SE,NW
const short eta[8] = {0, 0, 1, -1, 1, -1, 1, -1};
const short phi[8] = {1, -1, 0, 0, 1, -1, -1, 1};
for (size_t k = 0; k < 8; k++)
if (alpakaPfRecHits[i].neighbours()(k) != -1)
pfrh.addNeighbour(eta[k], phi[k], 0, alpakaPfRecHits[i].neighbours()(k));
}
}

event.emplace(legacyPfRecHitsToken_, out);
Expand Down
Original file line number Diff line number Diff line change
Expand Up @@ -39,10 +39,13 @@ namespace ALPAKA_ACCELERATOR_NAMESPACE {

reco::PFRecHitDeviceCollection pfRecHits{(int)num_recHits, event.queue()};

PFRecHitProducerKernel<CAL> kernel{event.queue(), num_recHits};
for (const auto& token : recHitsToken_)
kernel.processRecHits(event.queue(), event.get(token.first), setup.getData(token.second), topology, pfRecHits);
kernel.associateTopologyInfo(event.queue(), topology, pfRecHits);
if (num_recHits != 0) {
PFRecHitProducerKernel<CAL> kernel{event.queue(), num_recHits};
for (const auto& token : recHitsToken_)
kernel.processRecHits(
event.queue(), event.get(token.first), setup.getData(token.second), topology, pfRecHits);
kernel.associateTopologyInfo(event.queue(), topology, pfRecHits);
}

if (synchronise_)
alpaka::wait(event.queue());
Expand Down

0 comments on commit fc8422f

Please sign in to comment.