Skip to content
Snippets Groups Projects
Commit a59b6729 authored by Frank Winklmeier's avatar Frank Winklmeier
Browse files

Merge branch 'hepmc3_nightly_fixes_16122020_part_1' into 'master'

HepMC3 nightlies branch fixes 16122020 part 1

See merge request atlas/athena!39250
parents 5150ef06 dccda59d
No related branches found
No related tags found
No related merge requests found
......@@ -55,6 +55,7 @@ inline void line(std::ostream& os,const GenVertex& v) {v.print(os);}
inline void line(std::ostream& os,const GenVertex* v) {v->print(os);}
}
inline int barcode(ConstGenVertexPtr p) { return p->barcode();}
inline int barcode(const GenVertex p) { return p.barcode();}
inline void* raw_pointer(GenVertexPtr p) { return p;}
inline const void* raw_pointer(ConstGenVertexPtr p) { return p;}
inline std::ostream& operator<<( std::ostream& os, const GenVertex* v ) { if (v) return os<<(*v); else return os;}
......
......@@ -249,11 +249,16 @@ StatusCode TruthHitAnalysis::execute() {
if (currentGenEventIter != mcCollection->end()) {
int nvtx = 0;
int nvtx_sec=0;
for (HepMC::GenEvent::vertex_const_iterator vtx=(*currentGenEventIter)->vertices_begin(); vtx!=(*currentGenEventIter)->vertices_end(); ++vtx) {
double x = (*vtx)->position().x();
double y = (*vtx)->position().y();
double z = (*vtx)->position().z();
double r = sqrt(x*x+y*y);
#ifdef HEPMC3
for (auto vtx: (*currentGenEventIter)->vertices()) {
#else
for (HepMC::GenEvent::vertex_const_iterator vtxit=(*currentGenEventIter)->vertices_begin(); vtxit!=(*currentGenEventIter)->vertices_end(); ++vtxit) {
auto vtx=*vtxit;
#endif
double x = vtx->position().x();
double y = vtx->position().y();
double z = vtx->position().z();
double r = std::sqrt(x*x+y*y);
m_h_vtx_x->Fill(x);
m_h_vtx_y->Fill(y);
m_h_vtx_r->Fill(r);
......@@ -284,9 +289,9 @@ StatusCode TruthHitAnalysis::execute() {
int npart_prim=0;
int npart_sec=0;
HepMC::GenEvent::particle_const_iterator currentGenParticleIter;
for (currentGenParticleIter=(*currentGenEventIter)->particles_begin(); currentGenParticleIter!=(*currentGenEventIter)->particles_end(); ++currentGenParticleIter) {
const HepMC::FourVector mom = (*currentGenParticleIter)->momentum();
for (auto currentGenParticle: *(*currentGenEventIter)) {
const HepMC::FourVector mom = currentGenParticle->momentum();
m_h_truth_px->Fill(mom.x());
m_h_truth_py->Fill(mom.y());
......@@ -294,27 +299,27 @@ StatusCode TruthHitAnalysis::execute() {
m_h_truth_pt->Fill(mom.perp());
m_h_truth_eta->Fill(mom.eta());
m_h_truth_phi->Fill(mom.phi());
m_h_barcode->Fill((*currentGenParticleIter)->barcode());
m_h_part_status->Fill((*currentGenParticleIter)->status());
m_h_barcode->Fill(HepMC::barcode(currentGenParticle));
m_h_part_status->Fill(currentGenParticle->status());
m_truth_px->push_back(mom.x());
m_truth_py->push_back(mom.y());
m_truth_pz->push_back(mom.z());
m_truth_pt->push_back(mom.perp());
m_truth_eta->push_back(mom.eta());
m_truth_phi->push_back(mom.phi());
m_barcode->push_back((*currentGenParticleIter)->barcode());
m_status->push_back((*currentGenParticleIter)->status());
m_barcode->push_back(HepMC::barcode(currentGenParticle));
m_status->push_back(currentGenParticle->status());
int pdg = (*currentGenParticleIter)->pdg_id();
int pdg = currentGenParticle->pdg_id();
m_pdgid->push_back(pdg);
if ((*currentGenParticleIter)->barcode() < 200000) {
if (HepMC::barcode(currentGenParticle) < 200000) {
m_h_part_pdgid->Fill(pdg);
m_h_part_p->Fill(mom.rho());
m_h_part_p->Fill(std::sqrt(mom.x()*mom.x()+mom.y()*mom.y()+mom.z()*mom.z()));
m_h_part_eta->Fill(mom.eta());
m_h_part_phi->Fill(mom.phi());
++npart_prim;
if ((*currentGenParticleIter)->barcode() < 10000) {
if (HepMC::barcode(currentGenParticle) < 10000) {
m_h_n_generations->Fill(0);
}
else {
......@@ -324,7 +329,7 @@ StatusCode TruthHitAnalysis::execute() {
else {
m_h_part_pdgid_sec->Fill(pdg);
++npart_sec;
const int gen = (*currentGenParticleIter)->barcode()/1000000 + 2;
const int gen = HepMC::barcode(currentGenParticle)/1000000 + 2;
m_h_n_generations->Fill(gen);
}
} // End iteration over particles
......
......@@ -36,17 +36,17 @@ StatusCode DecayInFlyTruthTrajectoryBuilder::initialize() {
//================================================================
void DecayInFlyTruthTrajectoryBuilder::
buildTruthTrajectory(TruthTrajectory *result, const HepMC::GenParticle *input) const
buildTruthTrajectory(TruthTrajectory *result, HepMC::ConstGenParticlePtr input) const
{
result->clear();
if(input) {
const HepMC::GenParticle *next(nullptr);
const HepMC::GenParticle *current = input;
HepMC::ConstGenParticlePtr next{nullptr};
HepMC::ConstGenParticlePtr current = input;
// Extend trajectory outwards. The last particle should go at [0]
// in the TruthTrajectory, so we need to use a tmp storage while
// traversing the structure.
std::stack<const HepMC::GenParticle*> tmp;
std::stack<HepMC::ConstGenParticlePtr> tmp;
while( (next = getDaughter(current)) ) {
tmp.push(current = next);
}
......@@ -69,15 +69,20 @@ buildTruthTrajectory(TruthTrajectory *result, const HepMC::GenParticle *input) c
//================================================================
DecayInFlyTruthTrajectoryBuilder::MotherDaughter
DecayInFlyTruthTrajectoryBuilder::truthTrajectoryCuts(const HepMC::GenVertex *vtx) const
DecayInFlyTruthTrajectoryBuilder::truthTrajectoryCuts(HepMC::ConstGenVertexPtr vtx) const
{
const HepMC::GenParticle *mother(nullptr);
const HepMC::GenParticle *daughter(nullptr);
HepMC::GenParticlePtr mother{nullptr};
HepMC::GenParticlePtr daughter{nullptr};
// only truth vertices with 1 incoming particle
#ifdef HEPMC3
if(vtx && (vtx->particles_in().size() == 1)) {
mother = vtx->particles_in().front();
#else
if(vtx && (vtx->particles_in_size() == 1)) {
mother = *vtx->particles_in_const_begin();
#endif
// Allow status code 1 and 2. E.g. a pion that produced a long track can decay outside of InDet and have status==2.
if( mother && (mother->status() < 3) ) {
......@@ -92,11 +97,8 @@ DecayInFlyTruthTrajectoryBuilder::truthTrajectoryCuts(const HepMC::GenVertex *vt
if (vtx->particles_out_size() <= 2) {
int num_passed_cuts = 0;
const HepMC::GenParticle *passed_cuts(nullptr);
for(HepMC::GenVertex::particles_in_const_iterator it = vtx->particles_out_const_begin();
it != vtx->particles_out_const_end(); ++it) {
const HepMC::GenParticle *candidate = *it;
HepMC::GenParticlePtr passed_cuts{nullptr};
for(HepMC::GenParticlePtr candidate: *vtx){
if(candidate->pdg_id() == mother->pdg_id()) {
if(passed_cuts && (mother->pdg_id() == 11)) { // second negative electron is a special case
......@@ -128,9 +130,9 @@ DecayInFlyTruthTrajectoryBuilder::truthTrajectoryCuts(const HepMC::GenVertex *vt
}
//================================================================
const HepMC::GenParticle* DecayInFlyTruthTrajectoryBuilder::getDaughter(const HepMC::GenParticle* mother) const {
HepMC::ConstGenParticlePtr DecayInFlyTruthTrajectoryBuilder::getDaughter(HepMC::ConstGenParticlePtr mother) const {
const HepMC::GenParticle *daughter = nullptr;
HepMC::ConstGenParticlePtr daughter{nullptr};
if(mother) {
......@@ -145,9 +147,9 @@ const HepMC::GenParticle* DecayInFlyTruthTrajectoryBuilder::getDaughter(const He
}
//================================================================
const HepMC::GenParticle* DecayInFlyTruthTrajectoryBuilder::getMother(const HepMC::GenParticle* daughter) const {
HepMC::ConstGenParticlePtr DecayInFlyTruthTrajectoryBuilder::getMother(HepMC::ConstGenParticlePtr daughter) const {
const HepMC::GenParticle *mother = nullptr;
HepMC::ConstGenParticlePtr mother{nullptr};
if(daughter) {
......
......@@ -36,17 +36,17 @@ StatusCode ElasticTruthTrajectoryBuilder::initialize() {
//================================================================
void ElasticTruthTrajectoryBuilder::
buildTruthTrajectory(TruthTrajectory *result, const HepMC::GenParticle *input) const
buildTruthTrajectory(TruthTrajectory *result, HepMC::ConstGenParticlePtr input) const
{
result->clear();
if(input) {
const HepMC::GenParticle *next(nullptr);
const HepMC::GenParticle *current = input;
HepMC::ConstGenParticlePtr next(nullptr);
HepMC::ConstGenParticlePtr current = input;
// Extend trajectory outwards. The last particle should go at [0]
// in the TruthTrajectory, so we need to use a tmp storage while
// traversing the structure.
std::stack<const HepMC::GenParticle*> tmp;
std::stack<HepMC::ConstGenParticlePtr> tmp;
while( (next = getDaughter(current)) ) {
tmp.push(current = next);
}
......@@ -69,15 +69,19 @@ buildTruthTrajectory(TruthTrajectory *result, const HepMC::GenParticle *input) c
//================================================================
ElasticTruthTrajectoryBuilder::MotherDaughter
ElasticTruthTrajectoryBuilder::truthTrajectoryCuts(const HepMC::GenVertex *vtx) const
ElasticTruthTrajectoryBuilder::truthTrajectoryCuts(HepMC::ConstGenVertexPtr vtx) const
{
const HepMC::GenParticle *mother(nullptr);
const HepMC::GenParticle *daughter(nullptr);
HepMC::ConstGenParticlePtr mother{nullptr};
HepMC::ConstGenParticlePtr daughter{nullptr};
// only truth vertices with 1 incoming particle
#ifdef HEPMC3
if(vtx && (vtx->particles_in().size() == 1)) {
mother = vtx->particles_in().front();
#else
if(vtx && (vtx->particles_in_size() == 1)) {
mother = *vtx->particles_in_const_begin();
#endif
// Allow status code 1 and 2. E.g. a pion that produced a long track can decay outside of InDet and have status==2.
if( mother && (mother->status() < 3) ) {
......@@ -89,14 +93,14 @@ ElasticTruthTrajectoryBuilder::MotherDaughter
// is that with the higher energy (NOT pt).
//
// allow 1 outgoing to cover possible vertexes from interaction in detector material
#ifdef HEPMC3
if (vtx->particles_out().size() <= 2) {
#else
if (vtx->particles_out_size() <= 2) {
#endif
int num_passed_cuts = 0;
const HepMC::GenParticle *passed_cuts(nullptr);
for(HepMC::GenVertex::particles_in_const_iterator it = vtx->particles_out_const_begin();
it != vtx->particles_out_const_end(); ++it) {
const HepMC::GenParticle *candidate = *it;
HepMC::ConstGenParticlePtr passed_cuts{nullptr};
for(auto candidate: *vtx) {
if(candidate->pdg_id() == mother->pdg_id()) {
if(passed_cuts && (mother->pdg_id() == 11)) { // second negative electron is a special case
......@@ -122,9 +126,9 @@ ElasticTruthTrajectoryBuilder::MotherDaughter
}
//================================================================
const HepMC::GenParticle* ElasticTruthTrajectoryBuilder::getDaughter(const HepMC::GenParticle* mother) const {
HepMC::ConstGenParticlePtr ElasticTruthTrajectoryBuilder::getDaughter(HepMC::ConstGenParticlePtr mother) const {
const HepMC::GenParticle *daughter = nullptr;
HepMC::ConstGenParticlePtr daughter{nullptr};
if(mother) {
......@@ -139,9 +143,9 @@ const HepMC::GenParticle* ElasticTruthTrajectoryBuilder::getDaughter(const HepMC
}
//================================================================
const HepMC::GenParticle* ElasticTruthTrajectoryBuilder::getMother(const HepMC::GenParticle* daughter) const {
HepMC::ConstGenParticlePtr ElasticTruthTrajectoryBuilder::getMother(HepMC::ConstGenParticlePtr daughter) const {
const HepMC::GenParticle *mother = nullptr;
HepMC::ConstGenParticlePtr mother{nullptr};
if(daughter) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment