Commit 341f1f6a authored by Adam Edward Barton's avatar Adam Edward Barton
Browse files

Merge branch 'hepmc3_25.10.2021' into 'master'

Try to fix AGENE-2090

Closes AGENE-2090

See merge request atlas/athena!47586
parents fd9678ed 7f363146
/*
Copyright (C) 2002-2017 CERN for the benefit of the ATLAS collaboration
Copyright (C) 2002-2021 CERN for the benefit of the ATLAS collaboration
*/
#ifndef XAOD_ANALYSIS
......@@ -46,7 +46,6 @@ StatusCode FixHepMC::execute() {
// Event particle content cleaning -- remove "bad" structures
std::vector<HepMC::GenParticlePtr> toremove;
long seenThisEvent = 0;
/// @todo Use nicer particles accessor from TruthUtils / HepMC3 when it exists
for (auto ip: evt->particles()) {
// Skip this particle if (somehow) its pointer is null
if (!ip) continue;
......@@ -77,9 +76,30 @@ StatusCode FixHepMC::execute() {
ATH_MSG_DEBUG( "Found a bad particle in a decay chain : " );
if ( msgLvl( MSG::DEBUG ) ) HepMC::Print::line(ip);
}
/// Two types of bad particles: those w/o prod.vertex and with status other than 4.
if ( (!ip->production_vertex() || !ip->production_vertex()->id()) && ip->end_vertex() && ip->status() != 4 ) bad_particle = true;
/// Those w/o end vertex, but with bad status
if ( ip->production_vertex() && !ip->end_vertex() && ip->status() != 1 ) bad_particle = true;
// Only add to the toremove vector once, even if multiple tests match
if (bad_particle) toremove.push_back(ip);
}
// Some heuristics
// In case we have 3 particles, we try to add a vertex that correspond to 1->2 and 1->1 splitting.
if (toremove.size() == 3 || toremove.size() == 2) {
int no_endv = 0;
int no_prov = 0;
HepMC::FourVector sum(0,0,0,0);
for (auto part: toremove) if (!part->production_vertex() || !part->production_vertex()->id()) { no_prov++; sum += part->momentum();}
for (auto part: toremove) if (!part->end_vertex()) { no_endv++; sum -= part->momentum(); }
ATH_MSG_INFO("Heuristics: found " << toremove.size() << " bad particles. Try " << no_endv << "->" << no_prov << " splitting. The momenta sum is " << sum);
if (no_endv == 1 && (no_prov == 2 || no_prov == 1) && std::abs(sum.px()) < 1e-2 && std::abs(sum.py()) < 1e-2 && std::abs(sum.pz()) < 1e-2 ) {
auto v = HepMC::newGenVertexPtr();
for (auto part: toremove) if (!part->production_vertex() || !part->production_vertex()->id()) v->add_particle_out(part);
for (auto part: toremove) if (!part->end_vertex()) v->add_particle_in(part);
evt->add_vertex(v);
toremove.clear();
}
}
// Escape here if there's nothing more to do, otherwise do the cleaning
if (toremove.empty()) continue;
ATH_MSG_DEBUG("Cleaning event record of " << toremove.size() << " bad particles");
......@@ -112,7 +132,6 @@ StatusCode FixHepMC::execute() {
// Event particle content cleaning -- remove "bad" structures
std::vector<HepMC::GenParticlePtr> toremove; toremove.reserve(10);
long seenThisEvent = 0;
/// @todo Use nicer particles accessor from TruthUtils / HepMC3 when it exists
for (HepMC::GenEvent::particle_const_iterator ip = evt->particles_begin(); ip != evt->particles_end(); ++ip) {
// Skip this particle if (somehow) its pointer is null
if (*ip == NULL) continue;
......
......@@ -288,7 +288,8 @@ StatusCode TestHepMC::execute() {
beams.first=beams_t.at(0);
beams.second=beams_t.at(1);
} else {
ATH_MSG_WARNING("Invalid number of beam particles -- this generator interface should be fixed");
ATH_MSG_WARNING("Invalid number of beam particles " << beams_t.size() << " this generator interface should be fixed");
HepMC3::Print::content(*evt);
}
#else
auto beams = evt->beam_particles();
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment