ATLAS Offline Software
Loading...
Searching...
No Matches
MC::HepMC Namespace Reference

Namespaces

namespace  Print

Typedefs

typedef HepMC::GenVertex * GenVertexPtr
typedef const HepMC::GenVertex * ConstGenVertexPtr

Functions

template<class T>
int barcode (const T *p)
int barcode (int p)
template<class T, std::enable_if_t<!std::is_pointer< T >::value &&!std::is_same< T, int >::value, bool > = true>
int barcode (const T &p)
GenVertex::particles_out_const_iterator begin (const HepMC::GenVertex &v)
GenVertex::particles_out_const_iterator end (const HepMC::GenVertex &v)
GenVertexPtr newGenVertexPtr (const HepMC::FourVector &pos=HepMC::FourVector(0.0, 0.0, 0.0, 0.0), const int i=0)
int barcode_or_id (const ConstGenVertexPtr &p)
std::ostream & operator<< (std::ostream &os, const GenVertex *v)
bool set_ll_event_number (HepMC::GenEvent *e, long long int num)
long long int get_ll_event_number (const HepMC::GenEvent *e)
GenEvent::particle_iterator begin (HepMC::GenEvent &e)
GenEvent::particle_iterator end (HepMC::GenEvent &e)
GenEvent::particle_const_iterator begin (const HepMC::GenEvent &e)
GenEvent::particle_const_iterator end (const HepMC::GenEvent &e)
GenEvent * newGenEvent (const int a, const int b)
GenVertex * signal_process_vertex (const GenEvent *e)
void fillBarcodesAttribute (GenEvent *)
GenVertex * barcode_to_vertex (const GenEvent *e, int id)
GenParticlebarcode_to_particle (const GenEvent *e, int id)
int mpi (const GenEvent &e)
int mpi (const GenEvent *e)
int signal_process_id (const GenEvent &e)
int signal_process_id (const GenEvent *e)
void set_signal_process_id (GenEvent *e, const int i)
void set_mpi (GenEvent *e, const int i)
template<class T>
void set_random_states (GenEvent *e, std::vector< T > a)
template<class T>
void set_signal_process_vertex (GenEvent *e, T v)
GenEvent * copyemptyGenEvent (const GenEvent *inEvt)
template<class T>
bool suggest_barcode (T &p, int i)
template<class T>
bool suggest_barcode (T *p, int i)
template<>
bool suggest_barcode< std::unique_ptr< HepMC::GenParticle > > (std::unique_ptr< HepMC::GenParticle > &p, int i)
bool valid_beam_particles (const GenEvent *e)

Typedef Documentation

◆ ConstGenVertexPtr

typedef const HepMC::GenVertex* MC::HepMC::ConstGenVertexPtr

Definition at line 60 of file HepMCHelpers.h.

◆ GenVertexPtr

typedef HepMC::GenVertex* MC::HepMC::GenVertexPtr

Definition at line 59 of file HepMCHelpers.h.

Function Documentation

◆ barcode() [1/3]

template<class T, std::enable_if_t<!std::is_pointer< T >::value &&!std::is_same< T, int >::value, bool > = true>
int MC::HepMC::barcode ( const T & p)
inline

Definition at line 58 of file HepMCHelpers.h.

◆ barcode() [2/3]

template<class T>
int MC::HepMC::barcode ( const T * p)
inline

Definition at line 17 of file HepMCHelpers.h.

◆ barcode() [3/3]

int MC::HepMC::barcode ( int p)
inline

Definition at line 18 of file HepMCHelpers.h.

◆ barcode_or_id()

int MC::HepMC::barcode_or_id ( const ConstGenVertexPtr & p)
inline

Definition at line 71 of file HepMCHelpers.h.

◆ barcode_to_particle()

GenParticle * MC::HepMC::barcode_to_particle ( const GenEvent * e,
int id )
inline

Definition at line 628 of file HepMCHelpers.h.

649{
650inline
651auto particles_in (const HepMC::GenVertex* p) {
652 return std::ranges::subrange (p->particles_in_const_begin(),
653 p->particles_in_const_end());
654}
655}
656#endif
657
658namespace MC
659{
660 template <class VTX>
661 auto particles_in (const VTX* p) { return p->particles_in(); }
662 template <class VTX>
663 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
664
665 namespace Pythia8
666 {
668 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
669
670 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
671
672 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
673 }
674
675#include "AtlasPID.h"
676
678 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
679
681 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
682
684 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
685
687 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
688
690 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
691
693 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
694
696 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
697
699 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
700
702 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
703
705 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
706
710 template <class T> inline bool isStableOrSimDecayed(const T& p) {
711 const auto vertex = p->end_vertex();
712 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
713 }
714
716 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
717
719 template <class T> inline bool isSpecialNonInteracting(const T& p) {
720 const int apid = std::abs(p->pdg_id());
721 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
722 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
723 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
724 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
725 return false;
726 }
727
729
730 template <class T> T findMother(T thePart) {
731 auto partOriVert = thePart->production_vertex();
732 if (!partOriVert) return nullptr;
733
734 long partPDG = thePart->pdg_id();
735 long MotherPDG(0);
736
737 auto MothOriVert = partOriVert;
738 MothOriVert = nullptr;
739 T theMoth(nullptr);
740
741 size_t itr = 0;
742 do {
743 if (itr != 0) partOriVert = MothOriVert;
744 for ( const auto& p : particles_in(partOriVert) ) {
745 theMoth = p;
746 if (!theMoth) continue;
747 MotherPDG = theMoth->pdg_id();
748 MothOriVert = theMoth->production_vertex();
749 if (MotherPDG == partPDG) break;
750 }
751 itr++;
752 if (itr > 100) {
753 break;
754 }
755 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
756 MothOriVert != partOriVert);
757 return theMoth;
758 }
759
761
762 template <class C, class T> T findMatching(C TruthContainer, T p) {
763 T ptrPart = nullptr;
764 if (!p) return ptrPart;
765 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
766 for (T truthParticle : *TruthContainer) {
767 if (HepMC::is_sim_descendant(p,truthParticle)) {
768 ptrPart = truthParticle;
769 break;
770 }
771 }
772 }
773 else {
774 for (T truthParticle : TruthContainer) {
775 if (HepMC::is_sim_descendant(p,truthParticle)) {
776 ptrPart = truthParticle;
777 break;
778 }
779 }
780 }
781 return ptrPart;
782 }
784
785 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
786 auto prodVtx = thePart->production_vertex();
787 if (!prodVtx) return;
788 for (const auto& theMother: prodVtx->particles_in()) {
789 if (!theMother) continue;
790 allancestors.insert(theMother);
791 findParticleAncestors(theMother, allancestors);
792 }
793 }
794
796
797 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
798 auto endVtx = thePart->end_vertex();
799 if (!endVtx) return;
800 for (const auto& theDaughter: endVtx->particles_out()) {
801 if (!theDaughter) continue;
802 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
803 allstabledescendants.insert(theDaughter);
804 }
805 findParticleStableDescendants(theDaughter, allstabledescendants);
806 }
807 }
808
812
813 template <class T> bool isHardScatteringVertex(T pVert) {
814 if (pVert == nullptr) return false;
815 T pV = pVert;
816 int numOfPartIn(0);
817 int pdg(0);
818
819 do {
820 pVert = pV;
821 auto incoming = pVert->particles_in();
822 numOfPartIn = incoming.size();
823 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
824 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
825
826 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
827
828 if (numOfPartIn == 2) {
829 auto incoming = pVert->particles_in();
830 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
831 }
832 return false;
833}
834
838
839 template <class T, class U>
840 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
841 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
842 auto vtx = p->production_vertex();
843 if (!vtx) return false;
844 bool fromHad = false;
845 for ( const auto& parent : particles_in(vtx) ) {
846 if (!parent) continue;
847 // should this really go into parton-level territory?
848 // probably depends where BSM particles are being decayed
849 fromBSM |= isBSM(parent);
850 if (!isPhysical(parent)) return false;
851 fromTau |= isTau(parent);
852 if (isHadron(parent)&&!isBeam(parent)) {
853 if (!hadron) hadron = parent; // assumes linear hadron parentage
854 return true;
855 }
856 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
857 }
858 return fromHad;
859 }
860
863
864 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
865 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
866 decltype(thePart->end_vertex()) pVert(nullptr);
867 if (EndVert != nullptr) {
868 do {
869 bool samePart = false;
870 pVert = nullptr;
871 auto outgoing = EndVert->particles_out();
872 auto incoming = EndVert->particles_in();
873 for (const auto& itrDaug: outgoing) {
874 if (!itrDaug) continue;
875 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
876 // brem on generator level for tau
877 (outgoing.size() == 1 && incoming.size() == 1 &&
879 itrDaug->pdg_id() == thePart->pdg_id()) {
880 samePart = true;
881 pVert = itrDaug->end_vertex();
882 }
883 }
884 if (samePart) EndVert = pVert;
885 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
886 }
887 return EndVert;
888 }
889
891
892 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
893 if (!theVert) return {};
894 decltype(theVert->particles_out()) finalStatePart;
895 auto outgoing = theVert->particles_out();
896 for (const auto& thePart: outgoing) {
897 if (!thePart) continue;
898 finalStatePart.push_back(thePart);
899 if (isStable(thePart)) continue;
900 V pVert = findSimulatedEndVertex(thePart);
901 if (pVert == theVert) break; // to prevent Sherpa loop
902 if (pVert != nullptr) {
903 auto vecPart = findFinalStateParticles<V>(pVert);
904 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
905 }
906 }
907 return finalStatePart;
908 }
909#if !defined(XAOD_ANALYSIS)
910#include "AtlasHepMC/GenEvent.h"
911#ifdef HEPMC3
912inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
913inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
914#else
915inline void GeVToMeV(HepMC::GenEvent* evt) {
916 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
917 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
918 (*p)->momentum().py() * 1000,
919 (*p)->momentum().pz() * 1000,
920 (*p)->momentum().e() * 1000);
921 (*p)->set_momentum(fv);
922 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
923 }
924}
925inline void MeVToGeV(HepMC::GenEvent* evt) {
926 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
927 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
928 (*p)->momentum().py() / 1000,
929 (*p)->momentum().pz() / 1000,
930 (*p)->momentum().e() / 1000);
931 (*p)->set_momentum(fv);
932 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
933 }
934}
935#endif
936#endif
937}
938#endif
struct color C
bool is_same_generator_particle(const T1 &p1, const T2 &p2)
Method to establish if two particles in the GenEvent actually represent the same generated particle.
bool is_simulation_vertex(const T &v)
Method to establish if the vertex was created during simulation (TODO migrate to be based on status).
constexpr int SIM_STATUS_THRESHOLD
Constant definiting the status threshold for simulated particles, eg. can be used to separate generat...
int status(const T &p)
bool is_simulation_particle(const T &p)
Method to establish if a particle (or barcode) was created during the simulation (TODO update to be s...
constexpr bool is_smart_ptr_v
bool is_sim_descendant(const T1 &p1, const T2 &p2)
Method to check if the first particle is a descendant of the second in the simulation,...
bool isConditionB(const T &p)
bool isConditionA(const T &p)
To be understood.
bool isConditionC(const T &p)
T findMatching(C TruthContainer, T p)
Function to find a particle in container.
static const int GRAVITON
bool isZeroEnergyPhoton(const T &p)
Identify a photon with zero energy. Probably a workaround for a generator bug.
bool isHardScatteringVertex(T pVert)
Function to classify the vertex as hard scattering vertex.
void MeVToGeV(HepMC::GenEvent *evt)
bool isStableOrSimDecayed(const T &p)
Identify if particle is satble or decayed in simulation.
bool isPhoton(const T &p)
bool isSpecialNonInteracting(const T &p)
Identify a special non-interacting particles.
bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM)
Function to classify the particle.
bool isStable(const T &p)
Identify if the particle is stable, i.e. has not decayed.
bool isGeantino(const T &p)
bool isSimInteracting(const T &p)
Identify if the particle could interact with the detector during the simulation, e....
void GeVToMeV(HepMC::GenEvent *evt)
bool isEMInteracting(const T &p)
void findParticleAncestors(T thePart, std::set< T > &allancestors)
Function to find all ancestors of the particle.
bool isStrongInteracting(const T &p)
bool isInteracting(const T &p)
Identify if the particle with given PDG ID would not interact with the detector, i....
bool isMuon(const T &p)
bool isSUSY(const T &p)
auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out())
Function to find the stable particle descendants of the given vertex..
static const int NU_MU
bool isChargedNonShowering(const T &p)
Identify if the particle with given PDG ID would produce ID tracks but not shower in the detector if ...
bool isDecayed(const T &p)
Identify if the particle decayed.
T findMother(T thePart)
Function to get a mother of particle. MCTruthClassifier legacy.
auto particles_in(const HepMC::GenVertex *p)
void findParticleStableDescendants(T thePart, std::set< T > &allstabledescendants)
Function to get the particle stable MC daughters.
static const int NU_E
bool isBeam(const T &p)
Identify if the particle is beam particle.
static const int NU_TAU
bool isFinalState(const T &p)
Identify if the particle is final state particle.
bool isGenStable(const T &p)
Determine if the particle is stable at the generator (not det-sim) level,.
bool isHadron(const T &p)
bool isSimStable(const T &p)
Identify if the particle is considered stable at the post-detector-sim stage.
auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex())
Function to find the end vertex of a particle.
bool isTau(const T &p)
bool isPhysical(const T &p)
Identify if the particle is physical, i.e. is stable or decayed.
bool isBSM(const T &p)
APID: graviton and all Higgs extensions are BSM.
unsigned long long T

◆ barcode_to_vertex()

GenVertex * MC::HepMC::barcode_to_vertex ( const GenEvent * e,
int id )
inline

Definition at line 627 of file HepMCHelpers.h.

648{
649inline
650auto particles_in (const HepMC::GenVertex* p) {
651 return std::ranges::subrange (p->particles_in_const_begin(),
652 p->particles_in_const_end());
653}
654}
655#endif
656
657namespace MC
658{
659 template <class VTX>
660 auto particles_in (const VTX* p) { return p->particles_in(); }
661 template <class VTX>
662 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
663
664 namespace Pythia8
665 {
667 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
668
669 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
670
671 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
672 }
673
674#include "AtlasPID.h"
675
677 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
678
680 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
681
683 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
684
686 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
687
689 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
690
692 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
693
695 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
696
698 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
699
701 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
702
704 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
705
709 template <class T> inline bool isStableOrSimDecayed(const T& p) {
710 const auto vertex = p->end_vertex();
711 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
712 }
713
715 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
716
718 template <class T> inline bool isSpecialNonInteracting(const T& p) {
719 const int apid = std::abs(p->pdg_id());
720 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
721 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
722 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
723 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
724 return false;
725 }
726
728
729 template <class T> T findMother(T thePart) {
730 auto partOriVert = thePart->production_vertex();
731 if (!partOriVert) return nullptr;
732
733 long partPDG = thePart->pdg_id();
734 long MotherPDG(0);
735
736 auto MothOriVert = partOriVert;
737 MothOriVert = nullptr;
738 T theMoth(nullptr);
739
740 size_t itr = 0;
741 do {
742 if (itr != 0) partOriVert = MothOriVert;
743 for ( const auto& p : particles_in(partOriVert) ) {
744 theMoth = p;
745 if (!theMoth) continue;
746 MotherPDG = theMoth->pdg_id();
747 MothOriVert = theMoth->production_vertex();
748 if (MotherPDG == partPDG) break;
749 }
750 itr++;
751 if (itr > 100) {
752 break;
753 }
754 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
755 MothOriVert != partOriVert);
756 return theMoth;
757 }
758
760
761 template <class C, class T> T findMatching(C TruthContainer, T p) {
762 T ptrPart = nullptr;
763 if (!p) return ptrPart;
764 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
765 for (T truthParticle : *TruthContainer) {
766 if (HepMC::is_sim_descendant(p,truthParticle)) {
767 ptrPart = truthParticle;
768 break;
769 }
770 }
771 }
772 else {
773 for (T truthParticle : TruthContainer) {
774 if (HepMC::is_sim_descendant(p,truthParticle)) {
775 ptrPart = truthParticle;
776 break;
777 }
778 }
779 }
780 return ptrPart;
781 }
783
784 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
785 auto prodVtx = thePart->production_vertex();
786 if (!prodVtx) return;
787 for (const auto& theMother: prodVtx->particles_in()) {
788 if (!theMother) continue;
789 allancestors.insert(theMother);
790 findParticleAncestors(theMother, allancestors);
791 }
792 }
793
795
796 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
797 auto endVtx = thePart->end_vertex();
798 if (!endVtx) return;
799 for (const auto& theDaughter: endVtx->particles_out()) {
800 if (!theDaughter) continue;
801 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
802 allstabledescendants.insert(theDaughter);
803 }
804 findParticleStableDescendants(theDaughter, allstabledescendants);
805 }
806 }
807
811
812 template <class T> bool isHardScatteringVertex(T pVert) {
813 if (pVert == nullptr) return false;
814 T pV = pVert;
815 int numOfPartIn(0);
816 int pdg(0);
817
818 do {
819 pVert = pV;
820 auto incoming = pVert->particles_in();
821 numOfPartIn = incoming.size();
822 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
823 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
824
825 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
826
827 if (numOfPartIn == 2) {
828 auto incoming = pVert->particles_in();
829 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
830 }
831 return false;
832}
833
837
838 template <class T, class U>
839 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
840 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
841 auto vtx = p->production_vertex();
842 if (!vtx) return false;
843 bool fromHad = false;
844 for ( const auto& parent : particles_in(vtx) ) {
845 if (!parent) continue;
846 // should this really go into parton-level territory?
847 // probably depends where BSM particles are being decayed
848 fromBSM |= isBSM(parent);
849 if (!isPhysical(parent)) return false;
850 fromTau |= isTau(parent);
851 if (isHadron(parent)&&!isBeam(parent)) {
852 if (!hadron) hadron = parent; // assumes linear hadron parentage
853 return true;
854 }
855 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
856 }
857 return fromHad;
858 }
859
862
863 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
864 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
865 decltype(thePart->end_vertex()) pVert(nullptr);
866 if (EndVert != nullptr) {
867 do {
868 bool samePart = false;
869 pVert = nullptr;
870 auto outgoing = EndVert->particles_out();
871 auto incoming = EndVert->particles_in();
872 for (const auto& itrDaug: outgoing) {
873 if (!itrDaug) continue;
874 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
875 // brem on generator level for tau
876 (outgoing.size() == 1 && incoming.size() == 1 &&
878 itrDaug->pdg_id() == thePart->pdg_id()) {
879 samePart = true;
880 pVert = itrDaug->end_vertex();
881 }
882 }
883 if (samePart) EndVert = pVert;
884 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
885 }
886 return EndVert;
887 }
888
890
891 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
892 if (!theVert) return {};
893 decltype(theVert->particles_out()) finalStatePart;
894 auto outgoing = theVert->particles_out();
895 for (const auto& thePart: outgoing) {
896 if (!thePart) continue;
897 finalStatePart.push_back(thePart);
898 if (isStable(thePart)) continue;
899 V pVert = findSimulatedEndVertex(thePart);
900 if (pVert == theVert) break; // to prevent Sherpa loop
901 if (pVert != nullptr) {
902 auto vecPart = findFinalStateParticles<V>(pVert);
903 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
904 }
905 }
906 return finalStatePart;
907 }
908#if !defined(XAOD_ANALYSIS)
909#include "AtlasHepMC/GenEvent.h"
910#ifdef HEPMC3
911inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
912inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
913#else
914inline void GeVToMeV(HepMC::GenEvent* evt) {
915 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
916 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
917 (*p)->momentum().py() * 1000,
918 (*p)->momentum().pz() * 1000,
919 (*p)->momentum().e() * 1000);
920 (*p)->set_momentum(fv);
921 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
922 }
923}
924inline void MeVToGeV(HepMC::GenEvent* evt) {
925 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
926 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
927 (*p)->momentum().py() / 1000,
928 (*p)->momentum().pz() / 1000,
929 (*p)->momentum().e() / 1000);
930 (*p)->set_momentum(fv);
931 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
932 }
933}
934#endif
935#endif
936}
937#endif

◆ begin() [1/3]

GenEvent::particle_const_iterator MC::HepMC::begin ( const HepMC::GenEvent & e)
inline

Definition at line 622 of file HepMCHelpers.h.

643{
644inline
645auto particles_in (const HepMC::GenVertex* p) {
646 return std::ranges::subrange (p->particles_in_const_begin(),
647 p->particles_in_const_end());
648}
649}
650#endif
651
652namespace MC
653{
654 template <class VTX>
655 auto particles_in (const VTX* p) { return p->particles_in(); }
656 template <class VTX>
657 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
658
659 namespace Pythia8
660 {
662 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
663
664 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
665
666 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
667 }
668
669#include "AtlasPID.h"
670
672 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
673
675 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
676
678 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
679
681 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
682
684 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
685
687 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
688
690 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
691
693 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
694
696 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
697
699 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
700
704 template <class T> inline bool isStableOrSimDecayed(const T& p) {
705 const auto vertex = p->end_vertex();
706 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
707 }
708
710 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
711
713 template <class T> inline bool isSpecialNonInteracting(const T& p) {
714 const int apid = std::abs(p->pdg_id());
715 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
716 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
717 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
718 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
719 return false;
720 }
721
723
724 template <class T> T findMother(T thePart) {
725 auto partOriVert = thePart->production_vertex();
726 if (!partOriVert) return nullptr;
727
728 long partPDG = thePart->pdg_id();
729 long MotherPDG(0);
730
731 auto MothOriVert = partOriVert;
732 MothOriVert = nullptr;
733 T theMoth(nullptr);
734
735 size_t itr = 0;
736 do {
737 if (itr != 0) partOriVert = MothOriVert;
738 for ( const auto& p : particles_in(partOriVert) ) {
739 theMoth = p;
740 if (!theMoth) continue;
741 MotherPDG = theMoth->pdg_id();
742 MothOriVert = theMoth->production_vertex();
743 if (MotherPDG == partPDG) break;
744 }
745 itr++;
746 if (itr > 100) {
747 break;
748 }
749 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
750 MothOriVert != partOriVert);
751 return theMoth;
752 }
753
755
756 template <class C, class T> T findMatching(C TruthContainer, T p) {
757 T ptrPart = nullptr;
758 if (!p) return ptrPart;
759 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
760 for (T truthParticle : *TruthContainer) {
761 if (HepMC::is_sim_descendant(p,truthParticle)) {
762 ptrPart = truthParticle;
763 break;
764 }
765 }
766 }
767 else {
768 for (T truthParticle : TruthContainer) {
769 if (HepMC::is_sim_descendant(p,truthParticle)) {
770 ptrPart = truthParticle;
771 break;
772 }
773 }
774 }
775 return ptrPart;
776 }
778
779 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
780 auto prodVtx = thePart->production_vertex();
781 if (!prodVtx) return;
782 for (const auto& theMother: prodVtx->particles_in()) {
783 if (!theMother) continue;
784 allancestors.insert(theMother);
785 findParticleAncestors(theMother, allancestors);
786 }
787 }
788
790
791 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
792 auto endVtx = thePart->end_vertex();
793 if (!endVtx) return;
794 for (const auto& theDaughter: endVtx->particles_out()) {
795 if (!theDaughter) continue;
796 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
797 allstabledescendants.insert(theDaughter);
798 }
799 findParticleStableDescendants(theDaughter, allstabledescendants);
800 }
801 }
802
806
807 template <class T> bool isHardScatteringVertex(T pVert) {
808 if (pVert == nullptr) return false;
809 T pV = pVert;
810 int numOfPartIn(0);
811 int pdg(0);
812
813 do {
814 pVert = pV;
815 auto incoming = pVert->particles_in();
816 numOfPartIn = incoming.size();
817 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
818 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
819
820 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
821
822 if (numOfPartIn == 2) {
823 auto incoming = pVert->particles_in();
824 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
825 }
826 return false;
827}
828
832
833 template <class T, class U>
834 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
835 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
836 auto vtx = p->production_vertex();
837 if (!vtx) return false;
838 bool fromHad = false;
839 for ( const auto& parent : particles_in(vtx) ) {
840 if (!parent) continue;
841 // should this really go into parton-level territory?
842 // probably depends where BSM particles are being decayed
843 fromBSM |= isBSM(parent);
844 if (!isPhysical(parent)) return false;
845 fromTau |= isTau(parent);
846 if (isHadron(parent)&&!isBeam(parent)) {
847 if (!hadron) hadron = parent; // assumes linear hadron parentage
848 return true;
849 }
850 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
851 }
852 return fromHad;
853 }
854
857
858 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
859 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
860 decltype(thePart->end_vertex()) pVert(nullptr);
861 if (EndVert != nullptr) {
862 do {
863 bool samePart = false;
864 pVert = nullptr;
865 auto outgoing = EndVert->particles_out();
866 auto incoming = EndVert->particles_in();
867 for (const auto& itrDaug: outgoing) {
868 if (!itrDaug) continue;
869 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
870 // brem on generator level for tau
871 (outgoing.size() == 1 && incoming.size() == 1 &&
873 itrDaug->pdg_id() == thePart->pdg_id()) {
874 samePart = true;
875 pVert = itrDaug->end_vertex();
876 }
877 }
878 if (samePart) EndVert = pVert;
879 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
880 }
881 return EndVert;
882 }
883
885
886 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
887 if (!theVert) return {};
888 decltype(theVert->particles_out()) finalStatePart;
889 auto outgoing = theVert->particles_out();
890 for (const auto& thePart: outgoing) {
891 if (!thePart) continue;
892 finalStatePart.push_back(thePart);
893 if (isStable(thePart)) continue;
894 V pVert = findSimulatedEndVertex(thePart);
895 if (pVert == theVert) break; // to prevent Sherpa loop
896 if (pVert != nullptr) {
897 auto vecPart = findFinalStateParticles<V>(pVert);
898 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
899 }
900 }
901 return finalStatePart;
902 }
903#if !defined(XAOD_ANALYSIS)
904#include "AtlasHepMC/GenEvent.h"
905#ifdef HEPMC3
906inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
907inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
908#else
909inline void GeVToMeV(HepMC::GenEvent* evt) {
910 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
911 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
912 (*p)->momentum().py() * 1000,
913 (*p)->momentum().pz() * 1000,
914 (*p)->momentum().e() * 1000);
915 (*p)->set_momentum(fv);
916 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
917 }
918}
919inline void MeVToGeV(HepMC::GenEvent* evt) {
920 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
921 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
922 (*p)->momentum().py() / 1000,
923 (*p)->momentum().pz() / 1000,
924 (*p)->momentum().e() / 1000);
925 (*p)->set_momentum(fv);
926 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
927 }
928}
929#endif
930#endif
931}
932#endif

◆ begin() [2/3]

GenVertex::particles_out_const_iterator MC::HepMC::begin ( const HepMC::GenVertex & v)
inline

Definition at line 61 of file HepMCHelpers.h.

◆ begin() [3/3]

GenEvent::particle_iterator MC::HepMC::begin ( HepMC::GenEvent & e)
inline

Definition at line 620 of file HepMCHelpers.h.

641{
642inline
643auto particles_in (const HepMC::GenVertex* p) {
644 return std::ranges::subrange (p->particles_in_const_begin(),
645 p->particles_in_const_end());
646}
647}
648#endif
649
650namespace MC
651{
652 template <class VTX>
653 auto particles_in (const VTX* p) { return p->particles_in(); }
654 template <class VTX>
655 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
656
657 namespace Pythia8
658 {
660 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
661
662 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
663
664 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
665 }
666
667#include "AtlasPID.h"
668
670 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
671
673 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
674
676 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
677
679 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
680
682 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
683
685 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
686
688 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
689
691 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
692
694 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
695
697 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
698
702 template <class T> inline bool isStableOrSimDecayed(const T& p) {
703 const auto vertex = p->end_vertex();
704 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
705 }
706
708 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
709
711 template <class T> inline bool isSpecialNonInteracting(const T& p) {
712 const int apid = std::abs(p->pdg_id());
713 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
714 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
715 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
716 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
717 return false;
718 }
719
721
722 template <class T> T findMother(T thePart) {
723 auto partOriVert = thePart->production_vertex();
724 if (!partOriVert) return nullptr;
725
726 long partPDG = thePart->pdg_id();
727 long MotherPDG(0);
728
729 auto MothOriVert = partOriVert;
730 MothOriVert = nullptr;
731 T theMoth(nullptr);
732
733 size_t itr = 0;
734 do {
735 if (itr != 0) partOriVert = MothOriVert;
736 for ( const auto& p : particles_in(partOriVert) ) {
737 theMoth = p;
738 if (!theMoth) continue;
739 MotherPDG = theMoth->pdg_id();
740 MothOriVert = theMoth->production_vertex();
741 if (MotherPDG == partPDG) break;
742 }
743 itr++;
744 if (itr > 100) {
745 break;
746 }
747 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
748 MothOriVert != partOriVert);
749 return theMoth;
750 }
751
753
754 template <class C, class T> T findMatching(C TruthContainer, T p) {
755 T ptrPart = nullptr;
756 if (!p) return ptrPart;
757 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
758 for (T truthParticle : *TruthContainer) {
759 if (HepMC::is_sim_descendant(p,truthParticle)) {
760 ptrPart = truthParticle;
761 break;
762 }
763 }
764 }
765 else {
766 for (T truthParticle : TruthContainer) {
767 if (HepMC::is_sim_descendant(p,truthParticle)) {
768 ptrPart = truthParticle;
769 break;
770 }
771 }
772 }
773 return ptrPart;
774 }
776
777 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
778 auto prodVtx = thePart->production_vertex();
779 if (!prodVtx) return;
780 for (const auto& theMother: prodVtx->particles_in()) {
781 if (!theMother) continue;
782 allancestors.insert(theMother);
783 findParticleAncestors(theMother, allancestors);
784 }
785 }
786
788
789 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
790 auto endVtx = thePart->end_vertex();
791 if (!endVtx) return;
792 for (const auto& theDaughter: endVtx->particles_out()) {
793 if (!theDaughter) continue;
794 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
795 allstabledescendants.insert(theDaughter);
796 }
797 findParticleStableDescendants(theDaughter, allstabledescendants);
798 }
799 }
800
804
805 template <class T> bool isHardScatteringVertex(T pVert) {
806 if (pVert == nullptr) return false;
807 T pV = pVert;
808 int numOfPartIn(0);
809 int pdg(0);
810
811 do {
812 pVert = pV;
813 auto incoming = pVert->particles_in();
814 numOfPartIn = incoming.size();
815 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
816 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
817
818 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
819
820 if (numOfPartIn == 2) {
821 auto incoming = pVert->particles_in();
822 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
823 }
824 return false;
825}
826
830
831 template <class T, class U>
832 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
833 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
834 auto vtx = p->production_vertex();
835 if (!vtx) return false;
836 bool fromHad = false;
837 for ( const auto& parent : particles_in(vtx) ) {
838 if (!parent) continue;
839 // should this really go into parton-level territory?
840 // probably depends where BSM particles are being decayed
841 fromBSM |= isBSM(parent);
842 if (!isPhysical(parent)) return false;
843 fromTau |= isTau(parent);
844 if (isHadron(parent)&&!isBeam(parent)) {
845 if (!hadron) hadron = parent; // assumes linear hadron parentage
846 return true;
847 }
848 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
849 }
850 return fromHad;
851 }
852
855
856 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
857 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
858 decltype(thePart->end_vertex()) pVert(nullptr);
859 if (EndVert != nullptr) {
860 do {
861 bool samePart = false;
862 pVert = nullptr;
863 auto outgoing = EndVert->particles_out();
864 auto incoming = EndVert->particles_in();
865 for (const auto& itrDaug: outgoing) {
866 if (!itrDaug) continue;
867 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
868 // brem on generator level for tau
869 (outgoing.size() == 1 && incoming.size() == 1 &&
871 itrDaug->pdg_id() == thePart->pdg_id()) {
872 samePart = true;
873 pVert = itrDaug->end_vertex();
874 }
875 }
876 if (samePart) EndVert = pVert;
877 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
878 }
879 return EndVert;
880 }
881
883
884 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
885 if (!theVert) return {};
886 decltype(theVert->particles_out()) finalStatePart;
887 auto outgoing = theVert->particles_out();
888 for (const auto& thePart: outgoing) {
889 if (!thePart) continue;
890 finalStatePart.push_back(thePart);
891 if (isStable(thePart)) continue;
892 V pVert = findSimulatedEndVertex(thePart);
893 if (pVert == theVert) break; // to prevent Sherpa loop
894 if (pVert != nullptr) {
895 auto vecPart = findFinalStateParticles<V>(pVert);
896 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
897 }
898 }
899 return finalStatePart;
900 }
901#if !defined(XAOD_ANALYSIS)
902#include "AtlasHepMC/GenEvent.h"
903#ifdef HEPMC3
904inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
905inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
906#else
907inline void GeVToMeV(HepMC::GenEvent* evt) {
908 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
909 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
910 (*p)->momentum().py() * 1000,
911 (*p)->momentum().pz() * 1000,
912 (*p)->momentum().e() * 1000);
913 (*p)->set_momentum(fv);
914 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
915 }
916}
917inline void MeVToGeV(HepMC::GenEvent* evt) {
918 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
919 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
920 (*p)->momentum().py() / 1000,
921 (*p)->momentum().pz() / 1000,
922 (*p)->momentum().e() / 1000);
923 (*p)->set_momentum(fv);
924 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
925 }
926}
927#endif
928#endif
929}
930#endif

◆ copyemptyGenEvent()

GenEvent * MC::HepMC::copyemptyGenEvent ( const GenEvent * inEvt)
inline

Definition at line 653 of file HepMCHelpers.h.

674 {
675inline
676auto particles_in (const HepMC::GenVertex* p) {
677 return std::ranges::subrange (p->particles_in_const_begin(),
678 p->particles_in_const_end());
679}
680}
681#endif
682
683namespace MC
684{
685 template <class VTX>
686 auto particles_in (const VTX* p) { return p->particles_in(); }
687 template <class VTX>
688 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
689
690 namespace Pythia8
691 {
693 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
694
695 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
696
697 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
698 }
699
700#include "AtlasPID.h"
701
703 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
704
706 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
707
709 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
710
712 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
713
715 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
716
718 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
719
721 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
722
724 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
725
727 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
728
730 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
731
735 template <class T> inline bool isStableOrSimDecayed(const T& p) {
736 const auto vertex = p->end_vertex();
737 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
738 }
739
741 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
742
744 template <class T> inline bool isSpecialNonInteracting(const T& p) {
745 const int apid = std::abs(p->pdg_id());
746 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
747 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
748 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
749 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
750 return false;
751 }
752
754
755 template <class T> T findMother(T thePart) {
756 auto partOriVert = thePart->production_vertex();
757 if (!partOriVert) return nullptr;
758
759 long partPDG = thePart->pdg_id();
760 long MotherPDG(0);
761
762 auto MothOriVert = partOriVert;
763 MothOriVert = nullptr;
764 T theMoth(nullptr);
765
766 size_t itr = 0;
767 do {
768 if (itr != 0) partOriVert = MothOriVert;
769 for ( const auto& p : particles_in(partOriVert) ) {
770 theMoth = p;
771 if (!theMoth) continue;
772 MotherPDG = theMoth->pdg_id();
773 MothOriVert = theMoth->production_vertex();
774 if (MotherPDG == partPDG) break;
775 }
776 itr++;
777 if (itr > 100) {
778 break;
779 }
780 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
781 MothOriVert != partOriVert);
782 return theMoth;
783 }
784
786
787 template <class C, class T> T findMatching(C TruthContainer, T p) {
788 T ptrPart = nullptr;
789 if (!p) return ptrPart;
790 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
791 for (T truthParticle : *TruthContainer) {
792 if (HepMC::is_sim_descendant(p,truthParticle)) {
793 ptrPart = truthParticle;
794 break;
795 }
796 }
797 }
798 else {
799 for (T truthParticle : TruthContainer) {
800 if (HepMC::is_sim_descendant(p,truthParticle)) {
801 ptrPart = truthParticle;
802 break;
803 }
804 }
805 }
806 return ptrPart;
807 }
809
810 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
811 auto prodVtx = thePart->production_vertex();
812 if (!prodVtx) return;
813 for (const auto& theMother: prodVtx->particles_in()) {
814 if (!theMother) continue;
815 allancestors.insert(theMother);
816 findParticleAncestors(theMother, allancestors);
817 }
818 }
819
821
822 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
823 auto endVtx = thePart->end_vertex();
824 if (!endVtx) return;
825 for (const auto& theDaughter: endVtx->particles_out()) {
826 if (!theDaughter) continue;
827 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
828 allstabledescendants.insert(theDaughter);
829 }
830 findParticleStableDescendants(theDaughter, allstabledescendants);
831 }
832 }
833
837
838 template <class T> bool isHardScatteringVertex(T pVert) {
839 if (pVert == nullptr) return false;
840 T pV = pVert;
841 int numOfPartIn(0);
842 int pdg(0);
843
844 do {
845 pVert = pV;
846 auto incoming = pVert->particles_in();
847 numOfPartIn = incoming.size();
848 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
849 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
850
851 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
852
853 if (numOfPartIn == 2) {
854 auto incoming = pVert->particles_in();
855 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
856 }
857 return false;
858}
859
863
864 template <class T, class U>
865 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
866 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
867 auto vtx = p->production_vertex();
868 if (!vtx) return false;
869 bool fromHad = false;
870 for ( const auto& parent : particles_in(vtx) ) {
871 if (!parent) continue;
872 // should this really go into parton-level territory?
873 // probably depends where BSM particles are being decayed
874 fromBSM |= isBSM(parent);
875 if (!isPhysical(parent)) return false;
876 fromTau |= isTau(parent);
877 if (isHadron(parent)&&!isBeam(parent)) {
878 if (!hadron) hadron = parent; // assumes linear hadron parentage
879 return true;
880 }
881 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
882 }
883 return fromHad;
884 }
885
888
889 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
890 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
891 decltype(thePart->end_vertex()) pVert(nullptr);
892 if (EndVert != nullptr) {
893 do {
894 bool samePart = false;
895 pVert = nullptr;
896 auto outgoing = EndVert->particles_out();
897 auto incoming = EndVert->particles_in();
898 for (const auto& itrDaug: outgoing) {
899 if (!itrDaug) continue;
900 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
901 // brem on generator level for tau
902 (outgoing.size() == 1 && incoming.size() == 1 &&
904 itrDaug->pdg_id() == thePart->pdg_id()) {
905 samePart = true;
906 pVert = itrDaug->end_vertex();
907 }
908 }
909 if (samePart) EndVert = pVert;
910 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
911 }
912 return EndVert;
913 }
914
916
917 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
918 if (!theVert) return {};
919 decltype(theVert->particles_out()) finalStatePart;
920 auto outgoing = theVert->particles_out();
921 for (const auto& thePart: outgoing) {
922 if (!thePart) continue;
923 finalStatePart.push_back(thePart);
924 if (isStable(thePart)) continue;
925 V pVert = findSimulatedEndVertex(thePart);
926 if (pVert == theVert) break; // to prevent Sherpa loop
927 if (pVert != nullptr) {
928 auto vecPart = findFinalStateParticles<V>(pVert);
929 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
930 }
931 }
932 return finalStatePart;
933 }
934#if !defined(XAOD_ANALYSIS)
935#include "AtlasHepMC/GenEvent.h"
936#ifdef HEPMC3
937inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
938inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
939#else
940inline void GeVToMeV(HepMC::GenEvent* evt) {
941 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
942 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
943 (*p)->momentum().py() * 1000,
944 (*p)->momentum().pz() * 1000,
945 (*p)->momentum().e() * 1000);
946 (*p)->set_momentum(fv);
947 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
948 }
949}
950inline void MeVToGeV(HepMC::GenEvent* evt) {
951 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
952 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
953 (*p)->momentum().py() / 1000,
954 (*p)->momentum().pz() / 1000,
955 (*p)->momentum().e() / 1000);
956 (*p)->set_momentum(fv);
957 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
958 }
959}
960#endif
961#endif
962}
963#endif

◆ end() [1/3]

GenEvent::particle_const_iterator MC::HepMC::end ( const HepMC::GenEvent & e)
inline

Definition at line 623 of file HepMCHelpers.h.

644{
645inline
646auto particles_in (const HepMC::GenVertex* p) {
647 return std::ranges::subrange (p->particles_in_const_begin(),
648 p->particles_in_const_end());
649}
650}
651#endif
652
653namespace MC
654{
655 template <class VTX>
656 auto particles_in (const VTX* p) { return p->particles_in(); }
657 template <class VTX>
658 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
659
660 namespace Pythia8
661 {
663 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
664
665 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
666
667 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
668 }
669
670#include "AtlasPID.h"
671
673 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
674
676 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
677
679 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
680
682 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
683
685 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
686
688 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
689
691 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
692
694 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
695
697 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
698
700 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
701
705 template <class T> inline bool isStableOrSimDecayed(const T& p) {
706 const auto vertex = p->end_vertex();
707 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
708 }
709
711 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
712
714 template <class T> inline bool isSpecialNonInteracting(const T& p) {
715 const int apid = std::abs(p->pdg_id());
716 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
717 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
718 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
719 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
720 return false;
721 }
722
724
725 template <class T> T findMother(T thePart) {
726 auto partOriVert = thePart->production_vertex();
727 if (!partOriVert) return nullptr;
728
729 long partPDG = thePart->pdg_id();
730 long MotherPDG(0);
731
732 auto MothOriVert = partOriVert;
733 MothOriVert = nullptr;
734 T theMoth(nullptr);
735
736 size_t itr = 0;
737 do {
738 if (itr != 0) partOriVert = MothOriVert;
739 for ( const auto& p : particles_in(partOriVert) ) {
740 theMoth = p;
741 if (!theMoth) continue;
742 MotherPDG = theMoth->pdg_id();
743 MothOriVert = theMoth->production_vertex();
744 if (MotherPDG == partPDG) break;
745 }
746 itr++;
747 if (itr > 100) {
748 break;
749 }
750 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
751 MothOriVert != partOriVert);
752 return theMoth;
753 }
754
756
757 template <class C, class T> T findMatching(C TruthContainer, T p) {
758 T ptrPart = nullptr;
759 if (!p) return ptrPart;
760 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
761 for (T truthParticle : *TruthContainer) {
762 if (HepMC::is_sim_descendant(p,truthParticle)) {
763 ptrPart = truthParticle;
764 break;
765 }
766 }
767 }
768 else {
769 for (T truthParticle : TruthContainer) {
770 if (HepMC::is_sim_descendant(p,truthParticle)) {
771 ptrPart = truthParticle;
772 break;
773 }
774 }
775 }
776 return ptrPart;
777 }
779
780 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
781 auto prodVtx = thePart->production_vertex();
782 if (!prodVtx) return;
783 for (const auto& theMother: prodVtx->particles_in()) {
784 if (!theMother) continue;
785 allancestors.insert(theMother);
786 findParticleAncestors(theMother, allancestors);
787 }
788 }
789
791
792 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
793 auto endVtx = thePart->end_vertex();
794 if (!endVtx) return;
795 for (const auto& theDaughter: endVtx->particles_out()) {
796 if (!theDaughter) continue;
797 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
798 allstabledescendants.insert(theDaughter);
799 }
800 findParticleStableDescendants(theDaughter, allstabledescendants);
801 }
802 }
803
807
808 template <class T> bool isHardScatteringVertex(T pVert) {
809 if (pVert == nullptr) return false;
810 T pV = pVert;
811 int numOfPartIn(0);
812 int pdg(0);
813
814 do {
815 pVert = pV;
816 auto incoming = pVert->particles_in();
817 numOfPartIn = incoming.size();
818 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
819 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
820
821 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
822
823 if (numOfPartIn == 2) {
824 auto incoming = pVert->particles_in();
825 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
826 }
827 return false;
828}
829
833
834 template <class T, class U>
835 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
836 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
837 auto vtx = p->production_vertex();
838 if (!vtx) return false;
839 bool fromHad = false;
840 for ( const auto& parent : particles_in(vtx) ) {
841 if (!parent) continue;
842 // should this really go into parton-level territory?
843 // probably depends where BSM particles are being decayed
844 fromBSM |= isBSM(parent);
845 if (!isPhysical(parent)) return false;
846 fromTau |= isTau(parent);
847 if (isHadron(parent)&&!isBeam(parent)) {
848 if (!hadron) hadron = parent; // assumes linear hadron parentage
849 return true;
850 }
851 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
852 }
853 return fromHad;
854 }
855
858
859 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
860 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
861 decltype(thePart->end_vertex()) pVert(nullptr);
862 if (EndVert != nullptr) {
863 do {
864 bool samePart = false;
865 pVert = nullptr;
866 auto outgoing = EndVert->particles_out();
867 auto incoming = EndVert->particles_in();
868 for (const auto& itrDaug: outgoing) {
869 if (!itrDaug) continue;
870 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
871 // brem on generator level for tau
872 (outgoing.size() == 1 && incoming.size() == 1 &&
874 itrDaug->pdg_id() == thePart->pdg_id()) {
875 samePart = true;
876 pVert = itrDaug->end_vertex();
877 }
878 }
879 if (samePart) EndVert = pVert;
880 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
881 }
882 return EndVert;
883 }
884
886
887 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
888 if (!theVert) return {};
889 decltype(theVert->particles_out()) finalStatePart;
890 auto outgoing = theVert->particles_out();
891 for (const auto& thePart: outgoing) {
892 if (!thePart) continue;
893 finalStatePart.push_back(thePart);
894 if (isStable(thePart)) continue;
895 V pVert = findSimulatedEndVertex(thePart);
896 if (pVert == theVert) break; // to prevent Sherpa loop
897 if (pVert != nullptr) {
898 auto vecPart = findFinalStateParticles<V>(pVert);
899 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
900 }
901 }
902 return finalStatePart;
903 }
904#if !defined(XAOD_ANALYSIS)
905#include "AtlasHepMC/GenEvent.h"
906#ifdef HEPMC3
907inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
908inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
909#else
910inline void GeVToMeV(HepMC::GenEvent* evt) {
911 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
912 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
913 (*p)->momentum().py() * 1000,
914 (*p)->momentum().pz() * 1000,
915 (*p)->momentum().e() * 1000);
916 (*p)->set_momentum(fv);
917 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
918 }
919}
920inline void MeVToGeV(HepMC::GenEvent* evt) {
921 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
922 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
923 (*p)->momentum().py() / 1000,
924 (*p)->momentum().pz() / 1000,
925 (*p)->momentum().e() / 1000);
926 (*p)->set_momentum(fv);
927 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
928 }
929}
930#endif
931#endif
932}
933#endif

◆ end() [2/3]

GenVertex::particles_out_const_iterator MC::HepMC::end ( const HepMC::GenVertex & v)
inline

Definition at line 62 of file HepMCHelpers.h.

◆ end() [3/3]

GenEvent::particle_iterator MC::HepMC::end ( HepMC::GenEvent & e)
inline

Definition at line 621 of file HepMCHelpers.h.

642{
643inline
644auto particles_in (const HepMC::GenVertex* p) {
645 return std::ranges::subrange (p->particles_in_const_begin(),
646 p->particles_in_const_end());
647}
648}
649#endif
650
651namespace MC
652{
653 template <class VTX>
654 auto particles_in (const VTX* p) { return p->particles_in(); }
655 template <class VTX>
656 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
657
658 namespace Pythia8
659 {
661 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
662
663 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
664
665 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
666 }
667
668#include "AtlasPID.h"
669
671 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
672
674 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
675
677 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
678
680 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
681
683 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
684
686 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
687
689 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
690
692 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
693
695 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
696
698 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
699
703 template <class T> inline bool isStableOrSimDecayed(const T& p) {
704 const auto vertex = p->end_vertex();
705 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
706 }
707
709 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
710
712 template <class T> inline bool isSpecialNonInteracting(const T& p) {
713 const int apid = std::abs(p->pdg_id());
714 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
715 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
716 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
717 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
718 return false;
719 }
720
722
723 template <class T> T findMother(T thePart) {
724 auto partOriVert = thePart->production_vertex();
725 if (!partOriVert) return nullptr;
726
727 long partPDG = thePart->pdg_id();
728 long MotherPDG(0);
729
730 auto MothOriVert = partOriVert;
731 MothOriVert = nullptr;
732 T theMoth(nullptr);
733
734 size_t itr = 0;
735 do {
736 if (itr != 0) partOriVert = MothOriVert;
737 for ( const auto& p : particles_in(partOriVert) ) {
738 theMoth = p;
739 if (!theMoth) continue;
740 MotherPDG = theMoth->pdg_id();
741 MothOriVert = theMoth->production_vertex();
742 if (MotherPDG == partPDG) break;
743 }
744 itr++;
745 if (itr > 100) {
746 break;
747 }
748 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
749 MothOriVert != partOriVert);
750 return theMoth;
751 }
752
754
755 template <class C, class T> T findMatching(C TruthContainer, T p) {
756 T ptrPart = nullptr;
757 if (!p) return ptrPart;
758 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
759 for (T truthParticle : *TruthContainer) {
760 if (HepMC::is_sim_descendant(p,truthParticle)) {
761 ptrPart = truthParticle;
762 break;
763 }
764 }
765 }
766 else {
767 for (T truthParticle : TruthContainer) {
768 if (HepMC::is_sim_descendant(p,truthParticle)) {
769 ptrPart = truthParticle;
770 break;
771 }
772 }
773 }
774 return ptrPart;
775 }
777
778 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
779 auto prodVtx = thePart->production_vertex();
780 if (!prodVtx) return;
781 for (const auto& theMother: prodVtx->particles_in()) {
782 if (!theMother) continue;
783 allancestors.insert(theMother);
784 findParticleAncestors(theMother, allancestors);
785 }
786 }
787
789
790 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
791 auto endVtx = thePart->end_vertex();
792 if (!endVtx) return;
793 for (const auto& theDaughter: endVtx->particles_out()) {
794 if (!theDaughter) continue;
795 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
796 allstabledescendants.insert(theDaughter);
797 }
798 findParticleStableDescendants(theDaughter, allstabledescendants);
799 }
800 }
801
805
806 template <class T> bool isHardScatteringVertex(T pVert) {
807 if (pVert == nullptr) return false;
808 T pV = pVert;
809 int numOfPartIn(0);
810 int pdg(0);
811
812 do {
813 pVert = pV;
814 auto incoming = pVert->particles_in();
815 numOfPartIn = incoming.size();
816 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
817 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
818
819 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
820
821 if (numOfPartIn == 2) {
822 auto incoming = pVert->particles_in();
823 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
824 }
825 return false;
826}
827
831
832 template <class T, class U>
833 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
834 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
835 auto vtx = p->production_vertex();
836 if (!vtx) return false;
837 bool fromHad = false;
838 for ( const auto& parent : particles_in(vtx) ) {
839 if (!parent) continue;
840 // should this really go into parton-level territory?
841 // probably depends where BSM particles are being decayed
842 fromBSM |= isBSM(parent);
843 if (!isPhysical(parent)) return false;
844 fromTau |= isTau(parent);
845 if (isHadron(parent)&&!isBeam(parent)) {
846 if (!hadron) hadron = parent; // assumes linear hadron parentage
847 return true;
848 }
849 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
850 }
851 return fromHad;
852 }
853
856
857 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
858 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
859 decltype(thePart->end_vertex()) pVert(nullptr);
860 if (EndVert != nullptr) {
861 do {
862 bool samePart = false;
863 pVert = nullptr;
864 auto outgoing = EndVert->particles_out();
865 auto incoming = EndVert->particles_in();
866 for (const auto& itrDaug: outgoing) {
867 if (!itrDaug) continue;
868 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
869 // brem on generator level for tau
870 (outgoing.size() == 1 && incoming.size() == 1 &&
872 itrDaug->pdg_id() == thePart->pdg_id()) {
873 samePart = true;
874 pVert = itrDaug->end_vertex();
875 }
876 }
877 if (samePart) EndVert = pVert;
878 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
879 }
880 return EndVert;
881 }
882
884
885 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
886 if (!theVert) return {};
887 decltype(theVert->particles_out()) finalStatePart;
888 auto outgoing = theVert->particles_out();
889 for (const auto& thePart: outgoing) {
890 if (!thePart) continue;
891 finalStatePart.push_back(thePart);
892 if (isStable(thePart)) continue;
893 V pVert = findSimulatedEndVertex(thePart);
894 if (pVert == theVert) break; // to prevent Sherpa loop
895 if (pVert != nullptr) {
896 auto vecPart = findFinalStateParticles<V>(pVert);
897 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
898 }
899 }
900 return finalStatePart;
901 }
902#if !defined(XAOD_ANALYSIS)
903#include "AtlasHepMC/GenEvent.h"
904#ifdef HEPMC3
905inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
906inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
907#else
908inline void GeVToMeV(HepMC::GenEvent* evt) {
909 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
910 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
911 (*p)->momentum().py() * 1000,
912 (*p)->momentum().pz() * 1000,
913 (*p)->momentum().e() * 1000);
914 (*p)->set_momentum(fv);
915 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
916 }
917}
918inline void MeVToGeV(HepMC::GenEvent* evt) {
919 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
920 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
921 (*p)->momentum().py() / 1000,
922 (*p)->momentum().pz() / 1000,
923 (*p)->momentum().e() / 1000);
924 (*p)->set_momentum(fv);
925 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
926 }
927}
928#endif
929#endif
930}
931#endif

◆ fillBarcodesAttribute()

void MC::HepMC::fillBarcodesAttribute ( GenEvent * )
inline

Definition at line 626 of file HepMCHelpers.h.

647{
648inline
649auto particles_in (const HepMC::GenVertex* p) {
650 return std::ranges::subrange (p->particles_in_const_begin(),
651 p->particles_in_const_end());
652}
653}
654#endif
655
656namespace MC
657{
658 template <class VTX>
659 auto particles_in (const VTX* p) { return p->particles_in(); }
660 template <class VTX>
661 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
662
663 namespace Pythia8
664 {
666 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
667
668 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
669
670 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
671 }
672
673#include "AtlasPID.h"
674
676 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
677
679 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
680
682 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
683
685 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
686
688 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
689
691 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
692
694 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
695
697 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
698
700 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
701
703 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
704
708 template <class T> inline bool isStableOrSimDecayed(const T& p) {
709 const auto vertex = p->end_vertex();
710 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
711 }
712
714 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
715
717 template <class T> inline bool isSpecialNonInteracting(const T& p) {
718 const int apid = std::abs(p->pdg_id());
719 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
720 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
721 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
722 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
723 return false;
724 }
725
727
728 template <class T> T findMother(T thePart) {
729 auto partOriVert = thePart->production_vertex();
730 if (!partOriVert) return nullptr;
731
732 long partPDG = thePart->pdg_id();
733 long MotherPDG(0);
734
735 auto MothOriVert = partOriVert;
736 MothOriVert = nullptr;
737 T theMoth(nullptr);
738
739 size_t itr = 0;
740 do {
741 if (itr != 0) partOriVert = MothOriVert;
742 for ( const auto& p : particles_in(partOriVert) ) {
743 theMoth = p;
744 if (!theMoth) continue;
745 MotherPDG = theMoth->pdg_id();
746 MothOriVert = theMoth->production_vertex();
747 if (MotherPDG == partPDG) break;
748 }
749 itr++;
750 if (itr > 100) {
751 break;
752 }
753 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
754 MothOriVert != partOriVert);
755 return theMoth;
756 }
757
759
760 template <class C, class T> T findMatching(C TruthContainer, T p) {
761 T ptrPart = nullptr;
762 if (!p) return ptrPart;
763 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
764 for (T truthParticle : *TruthContainer) {
765 if (HepMC::is_sim_descendant(p,truthParticle)) {
766 ptrPart = truthParticle;
767 break;
768 }
769 }
770 }
771 else {
772 for (T truthParticle : TruthContainer) {
773 if (HepMC::is_sim_descendant(p,truthParticle)) {
774 ptrPart = truthParticle;
775 break;
776 }
777 }
778 }
779 return ptrPart;
780 }
782
783 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
784 auto prodVtx = thePart->production_vertex();
785 if (!prodVtx) return;
786 for (const auto& theMother: prodVtx->particles_in()) {
787 if (!theMother) continue;
788 allancestors.insert(theMother);
789 findParticleAncestors(theMother, allancestors);
790 }
791 }
792
794
795 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
796 auto endVtx = thePart->end_vertex();
797 if (!endVtx) return;
798 for (const auto& theDaughter: endVtx->particles_out()) {
799 if (!theDaughter) continue;
800 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
801 allstabledescendants.insert(theDaughter);
802 }
803 findParticleStableDescendants(theDaughter, allstabledescendants);
804 }
805 }
806
810
811 template <class T> bool isHardScatteringVertex(T pVert) {
812 if (pVert == nullptr) return false;
813 T pV = pVert;
814 int numOfPartIn(0);
815 int pdg(0);
816
817 do {
818 pVert = pV;
819 auto incoming = pVert->particles_in();
820 numOfPartIn = incoming.size();
821 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
822 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
823
824 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
825
826 if (numOfPartIn == 2) {
827 auto incoming = pVert->particles_in();
828 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
829 }
830 return false;
831}
832
836
837 template <class T, class U>
838 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
839 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
840 auto vtx = p->production_vertex();
841 if (!vtx) return false;
842 bool fromHad = false;
843 for ( const auto& parent : particles_in(vtx) ) {
844 if (!parent) continue;
845 // should this really go into parton-level territory?
846 // probably depends where BSM particles are being decayed
847 fromBSM |= isBSM(parent);
848 if (!isPhysical(parent)) return false;
849 fromTau |= isTau(parent);
850 if (isHadron(parent)&&!isBeam(parent)) {
851 if (!hadron) hadron = parent; // assumes linear hadron parentage
852 return true;
853 }
854 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
855 }
856 return fromHad;
857 }
858
861
862 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
863 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
864 decltype(thePart->end_vertex()) pVert(nullptr);
865 if (EndVert != nullptr) {
866 do {
867 bool samePart = false;
868 pVert = nullptr;
869 auto outgoing = EndVert->particles_out();
870 auto incoming = EndVert->particles_in();
871 for (const auto& itrDaug: outgoing) {
872 if (!itrDaug) continue;
873 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
874 // brem on generator level for tau
875 (outgoing.size() == 1 && incoming.size() == 1 &&
877 itrDaug->pdg_id() == thePart->pdg_id()) {
878 samePart = true;
879 pVert = itrDaug->end_vertex();
880 }
881 }
882 if (samePart) EndVert = pVert;
883 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
884 }
885 return EndVert;
886 }
887
889
890 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
891 if (!theVert) return {};
892 decltype(theVert->particles_out()) finalStatePart;
893 auto outgoing = theVert->particles_out();
894 for (const auto& thePart: outgoing) {
895 if (!thePart) continue;
896 finalStatePart.push_back(thePart);
897 if (isStable(thePart)) continue;
898 V pVert = findSimulatedEndVertex(thePart);
899 if (pVert == theVert) break; // to prevent Sherpa loop
900 if (pVert != nullptr) {
901 auto vecPart = findFinalStateParticles<V>(pVert);
902 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
903 }
904 }
905 return finalStatePart;
906 }
907#if !defined(XAOD_ANALYSIS)
908#include "AtlasHepMC/GenEvent.h"
909#ifdef HEPMC3
910inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
911inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
912#else
913inline void GeVToMeV(HepMC::GenEvent* evt) {
914 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
915 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
916 (*p)->momentum().py() * 1000,
917 (*p)->momentum().pz() * 1000,
918 (*p)->momentum().e() * 1000);
919 (*p)->set_momentum(fv);
920 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
921 }
922}
923inline void MeVToGeV(HepMC::GenEvent* evt) {
924 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
925 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
926 (*p)->momentum().py() / 1000,
927 (*p)->momentum().pz() / 1000,
928 (*p)->momentum().e() / 1000);
929 (*p)->set_momentum(fv);
930 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
931 }
932}
933#endif
934#endif
935}
936#endif

◆ get_ll_event_number()

long long int MC::HepMC::get_ll_event_number ( const HepMC::GenEvent * e)
inline

Definition at line 617 of file HepMCHelpers.h.

638 {
639inline
640auto particles_in (const HepMC::GenVertex* p) {
641 return std::ranges::subrange (p->particles_in_const_begin(),
642 p->particles_in_const_end());
643}
644}
645#endif
646
647namespace MC
648{
649 template <class VTX>
650 auto particles_in (const VTX* p) { return p->particles_in(); }
651 template <class VTX>
652 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
653
654 namespace Pythia8
655 {
657 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
658
659 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
660
661 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
662 }
663
664#include "AtlasPID.h"
665
667 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
668
670 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
671
673 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
674
676 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
677
679 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
680
682 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
683
685 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
686
688 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
689
691 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
692
694 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
695
699 template <class T> inline bool isStableOrSimDecayed(const T& p) {
700 const auto vertex = p->end_vertex();
701 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
702 }
703
705 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
706
708 template <class T> inline bool isSpecialNonInteracting(const T& p) {
709 const int apid = std::abs(p->pdg_id());
710 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
711 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
712 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
713 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
714 return false;
715 }
716
718
719 template <class T> T findMother(T thePart) {
720 auto partOriVert = thePart->production_vertex();
721 if (!partOriVert) return nullptr;
722
723 long partPDG = thePart->pdg_id();
724 long MotherPDG(0);
725
726 auto MothOriVert = partOriVert;
727 MothOriVert = nullptr;
728 T theMoth(nullptr);
729
730 size_t itr = 0;
731 do {
732 if (itr != 0) partOriVert = MothOriVert;
733 for ( const auto& p : particles_in(partOriVert) ) {
734 theMoth = p;
735 if (!theMoth) continue;
736 MotherPDG = theMoth->pdg_id();
737 MothOriVert = theMoth->production_vertex();
738 if (MotherPDG == partPDG) break;
739 }
740 itr++;
741 if (itr > 100) {
742 break;
743 }
744 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
745 MothOriVert != partOriVert);
746 return theMoth;
747 }
748
750
751 template <class C, class T> T findMatching(C TruthContainer, T p) {
752 T ptrPart = nullptr;
753 if (!p) return ptrPart;
754 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
755 for (T truthParticle : *TruthContainer) {
756 if (HepMC::is_sim_descendant(p,truthParticle)) {
757 ptrPart = truthParticle;
758 break;
759 }
760 }
761 }
762 else {
763 for (T truthParticle : TruthContainer) {
764 if (HepMC::is_sim_descendant(p,truthParticle)) {
765 ptrPart = truthParticle;
766 break;
767 }
768 }
769 }
770 return ptrPart;
771 }
773
774 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
775 auto prodVtx = thePart->production_vertex();
776 if (!prodVtx) return;
777 for (const auto& theMother: prodVtx->particles_in()) {
778 if (!theMother) continue;
779 allancestors.insert(theMother);
780 findParticleAncestors(theMother, allancestors);
781 }
782 }
783
785
786 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
787 auto endVtx = thePart->end_vertex();
788 if (!endVtx) return;
789 for (const auto& theDaughter: endVtx->particles_out()) {
790 if (!theDaughter) continue;
791 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
792 allstabledescendants.insert(theDaughter);
793 }
794 findParticleStableDescendants(theDaughter, allstabledescendants);
795 }
796 }
797
801
802 template <class T> bool isHardScatteringVertex(T pVert) {
803 if (pVert == nullptr) return false;
804 T pV = pVert;
805 int numOfPartIn(0);
806 int pdg(0);
807
808 do {
809 pVert = pV;
810 auto incoming = pVert->particles_in();
811 numOfPartIn = incoming.size();
812 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
813 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
814
815 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
816
817 if (numOfPartIn == 2) {
818 auto incoming = pVert->particles_in();
819 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
820 }
821 return false;
822}
823
827
828 template <class T, class U>
829 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
830 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
831 auto vtx = p->production_vertex();
832 if (!vtx) return false;
833 bool fromHad = false;
834 for ( const auto& parent : particles_in(vtx) ) {
835 if (!parent) continue;
836 // should this really go into parton-level territory?
837 // probably depends where BSM particles are being decayed
838 fromBSM |= isBSM(parent);
839 if (!isPhysical(parent)) return false;
840 fromTau |= isTau(parent);
841 if (isHadron(parent)&&!isBeam(parent)) {
842 if (!hadron) hadron = parent; // assumes linear hadron parentage
843 return true;
844 }
845 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
846 }
847 return fromHad;
848 }
849
852
853 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
854 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
855 decltype(thePart->end_vertex()) pVert(nullptr);
856 if (EndVert != nullptr) {
857 do {
858 bool samePart = false;
859 pVert = nullptr;
860 auto outgoing = EndVert->particles_out();
861 auto incoming = EndVert->particles_in();
862 for (const auto& itrDaug: outgoing) {
863 if (!itrDaug) continue;
864 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
865 // brem on generator level for tau
866 (outgoing.size() == 1 && incoming.size() == 1 &&
868 itrDaug->pdg_id() == thePart->pdg_id()) {
869 samePart = true;
870 pVert = itrDaug->end_vertex();
871 }
872 }
873 if (samePart) EndVert = pVert;
874 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
875 }
876 return EndVert;
877 }
878
880
881 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
882 if (!theVert) return {};
883 decltype(theVert->particles_out()) finalStatePart;
884 auto outgoing = theVert->particles_out();
885 for (const auto& thePart: outgoing) {
886 if (!thePart) continue;
887 finalStatePart.push_back(thePart);
888 if (isStable(thePart)) continue;
889 V pVert = findSimulatedEndVertex(thePart);
890 if (pVert == theVert) break; // to prevent Sherpa loop
891 if (pVert != nullptr) {
892 auto vecPart = findFinalStateParticles<V>(pVert);
893 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
894 }
895 }
896 return finalStatePart;
897 }
898#if !defined(XAOD_ANALYSIS)
899#include "AtlasHepMC/GenEvent.h"
900#ifdef HEPMC3
901inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
902inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
903#else
904inline void GeVToMeV(HepMC::GenEvent* evt) {
905 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
906 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
907 (*p)->momentum().py() * 1000,
908 (*p)->momentum().pz() * 1000,
909 (*p)->momentum().e() * 1000);
910 (*p)->set_momentum(fv);
911 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
912 }
913}
914inline void MeVToGeV(HepMC::GenEvent* evt) {
915 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
916 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
917 (*p)->momentum().py() / 1000,
918 (*p)->momentum().pz() / 1000,
919 (*p)->momentum().e() / 1000);
920 (*p)->set_momentum(fv);
921 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
922 }
923}
924#endif
925#endif
926}
927#endif

◆ mpi() [1/2]

int MC::HepMC::mpi ( const GenEvent & e)
inline

Definition at line 629 of file HepMCHelpers.h.

650 {
651inline
652auto particles_in (const HepMC::GenVertex* p) {
653 return std::ranges::subrange (p->particles_in_const_begin(),
654 p->particles_in_const_end());
655}
656}
657#endif
658
659namespace MC
660{
661 template <class VTX>
662 auto particles_in (const VTX* p) { return p->particles_in(); }
663 template <class VTX>
664 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
665
666 namespace Pythia8
667 {
669 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
670
671 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
672
673 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
674 }
675
676#include "AtlasPID.h"
677
679 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
680
682 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
683
685 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
686
688 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
689
691 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
692
694 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
695
697 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
698
700 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
701
703 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
704
706 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
707
711 template <class T> inline bool isStableOrSimDecayed(const T& p) {
712 const auto vertex = p->end_vertex();
713 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
714 }
715
717 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
718
720 template <class T> inline bool isSpecialNonInteracting(const T& p) {
721 const int apid = std::abs(p->pdg_id());
722 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
723 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
724 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
725 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
726 return false;
727 }
728
730
731 template <class T> T findMother(T thePart) {
732 auto partOriVert = thePart->production_vertex();
733 if (!partOriVert) return nullptr;
734
735 long partPDG = thePart->pdg_id();
736 long MotherPDG(0);
737
738 auto MothOriVert = partOriVert;
739 MothOriVert = nullptr;
740 T theMoth(nullptr);
741
742 size_t itr = 0;
743 do {
744 if (itr != 0) partOriVert = MothOriVert;
745 for ( const auto& p : particles_in(partOriVert) ) {
746 theMoth = p;
747 if (!theMoth) continue;
748 MotherPDG = theMoth->pdg_id();
749 MothOriVert = theMoth->production_vertex();
750 if (MotherPDG == partPDG) break;
751 }
752 itr++;
753 if (itr > 100) {
754 break;
755 }
756 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
757 MothOriVert != partOriVert);
758 return theMoth;
759 }
760
762
763 template <class C, class T> T findMatching(C TruthContainer, T p) {
764 T ptrPart = nullptr;
765 if (!p) return ptrPart;
766 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
767 for (T truthParticle : *TruthContainer) {
768 if (HepMC::is_sim_descendant(p,truthParticle)) {
769 ptrPart = truthParticle;
770 break;
771 }
772 }
773 }
774 else {
775 for (T truthParticle : TruthContainer) {
776 if (HepMC::is_sim_descendant(p,truthParticle)) {
777 ptrPart = truthParticle;
778 break;
779 }
780 }
781 }
782 return ptrPart;
783 }
785
786 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
787 auto prodVtx = thePart->production_vertex();
788 if (!prodVtx) return;
789 for (const auto& theMother: prodVtx->particles_in()) {
790 if (!theMother) continue;
791 allancestors.insert(theMother);
792 findParticleAncestors(theMother, allancestors);
793 }
794 }
795
797
798 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
799 auto endVtx = thePart->end_vertex();
800 if (!endVtx) return;
801 for (const auto& theDaughter: endVtx->particles_out()) {
802 if (!theDaughter) continue;
803 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
804 allstabledescendants.insert(theDaughter);
805 }
806 findParticleStableDescendants(theDaughter, allstabledescendants);
807 }
808 }
809
813
814 template <class T> bool isHardScatteringVertex(T pVert) {
815 if (pVert == nullptr) return false;
816 T pV = pVert;
817 int numOfPartIn(0);
818 int pdg(0);
819
820 do {
821 pVert = pV;
822 auto incoming = pVert->particles_in();
823 numOfPartIn = incoming.size();
824 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
825 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
826
827 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
828
829 if (numOfPartIn == 2) {
830 auto incoming = pVert->particles_in();
831 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
832 }
833 return false;
834}
835
839
840 template <class T, class U>
841 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
842 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
843 auto vtx = p->production_vertex();
844 if (!vtx) return false;
845 bool fromHad = false;
846 for ( const auto& parent : particles_in(vtx) ) {
847 if (!parent) continue;
848 // should this really go into parton-level territory?
849 // probably depends where BSM particles are being decayed
850 fromBSM |= isBSM(parent);
851 if (!isPhysical(parent)) return false;
852 fromTau |= isTau(parent);
853 if (isHadron(parent)&&!isBeam(parent)) {
854 if (!hadron) hadron = parent; // assumes linear hadron parentage
855 return true;
856 }
857 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
858 }
859 return fromHad;
860 }
861
864
865 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
866 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
867 decltype(thePart->end_vertex()) pVert(nullptr);
868 if (EndVert != nullptr) {
869 do {
870 bool samePart = false;
871 pVert = nullptr;
872 auto outgoing = EndVert->particles_out();
873 auto incoming = EndVert->particles_in();
874 for (const auto& itrDaug: outgoing) {
875 if (!itrDaug) continue;
876 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
877 // brem on generator level for tau
878 (outgoing.size() == 1 && incoming.size() == 1 &&
880 itrDaug->pdg_id() == thePart->pdg_id()) {
881 samePart = true;
882 pVert = itrDaug->end_vertex();
883 }
884 }
885 if (samePart) EndVert = pVert;
886 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
887 }
888 return EndVert;
889 }
890
892
893 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
894 if (!theVert) return {};
895 decltype(theVert->particles_out()) finalStatePart;
896 auto outgoing = theVert->particles_out();
897 for (const auto& thePart: outgoing) {
898 if (!thePart) continue;
899 finalStatePart.push_back(thePart);
900 if (isStable(thePart)) continue;
901 V pVert = findSimulatedEndVertex(thePart);
902 if (pVert == theVert) break; // to prevent Sherpa loop
903 if (pVert != nullptr) {
904 auto vecPart = findFinalStateParticles<V>(pVert);
905 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
906 }
907 }
908 return finalStatePart;
909 }
910#if !defined(XAOD_ANALYSIS)
911#include "AtlasHepMC/GenEvent.h"
912#ifdef HEPMC3
913inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
914inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
915#else
916inline void GeVToMeV(HepMC::GenEvent* evt) {
917 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
918 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
919 (*p)->momentum().py() * 1000,
920 (*p)->momentum().pz() * 1000,
921 (*p)->momentum().e() * 1000);
922 (*p)->set_momentum(fv);
923 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
924 }
925}
926inline void MeVToGeV(HepMC::GenEvent* evt) {
927 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
928 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
929 (*p)->momentum().py() / 1000,
930 (*p)->momentum().pz() / 1000,
931 (*p)->momentum().e() / 1000);
932 (*p)->set_momentum(fv);
933 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
934 }
935}
936#endif
937#endif
938}
939#endif

◆ mpi() [2/2]

int MC::HepMC::mpi ( const GenEvent * e)
inline

Definition at line 632 of file HepMCHelpers.h.

653 {
654inline
655auto particles_in (const HepMC::GenVertex* p) {
656 return std::ranges::subrange (p->particles_in_const_begin(),
657 p->particles_in_const_end());
658}
659}
660#endif
661
662namespace MC
663{
664 template <class VTX>
665 auto particles_in (const VTX* p) { return p->particles_in(); }
666 template <class VTX>
667 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
668
669 namespace Pythia8
670 {
672 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
673
674 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
675
676 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
677 }
678
679#include "AtlasPID.h"
680
682 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
683
685 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
686
688 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
689
691 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
692
694 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
695
697 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
698
700 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
701
703 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
704
706 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
707
709 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
710
714 template <class T> inline bool isStableOrSimDecayed(const T& p) {
715 const auto vertex = p->end_vertex();
716 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
717 }
718
720 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
721
723 template <class T> inline bool isSpecialNonInteracting(const T& p) {
724 const int apid = std::abs(p->pdg_id());
725 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
726 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
727 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
728 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
729 return false;
730 }
731
733
734 template <class T> T findMother(T thePart) {
735 auto partOriVert = thePart->production_vertex();
736 if (!partOriVert) return nullptr;
737
738 long partPDG = thePart->pdg_id();
739 long MotherPDG(0);
740
741 auto MothOriVert = partOriVert;
742 MothOriVert = nullptr;
743 T theMoth(nullptr);
744
745 size_t itr = 0;
746 do {
747 if (itr != 0) partOriVert = MothOriVert;
748 for ( const auto& p : particles_in(partOriVert) ) {
749 theMoth = p;
750 if (!theMoth) continue;
751 MotherPDG = theMoth->pdg_id();
752 MothOriVert = theMoth->production_vertex();
753 if (MotherPDG == partPDG) break;
754 }
755 itr++;
756 if (itr > 100) {
757 break;
758 }
759 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
760 MothOriVert != partOriVert);
761 return theMoth;
762 }
763
765
766 template <class C, class T> T findMatching(C TruthContainer, T p) {
767 T ptrPart = nullptr;
768 if (!p) return ptrPart;
769 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
770 for (T truthParticle : *TruthContainer) {
771 if (HepMC::is_sim_descendant(p,truthParticle)) {
772 ptrPart = truthParticle;
773 break;
774 }
775 }
776 }
777 else {
778 for (T truthParticle : TruthContainer) {
779 if (HepMC::is_sim_descendant(p,truthParticle)) {
780 ptrPart = truthParticle;
781 break;
782 }
783 }
784 }
785 return ptrPart;
786 }
788
789 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
790 auto prodVtx = thePart->production_vertex();
791 if (!prodVtx) return;
792 for (const auto& theMother: prodVtx->particles_in()) {
793 if (!theMother) continue;
794 allancestors.insert(theMother);
795 findParticleAncestors(theMother, allancestors);
796 }
797 }
798
800
801 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
802 auto endVtx = thePart->end_vertex();
803 if (!endVtx) return;
804 for (const auto& theDaughter: endVtx->particles_out()) {
805 if (!theDaughter) continue;
806 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
807 allstabledescendants.insert(theDaughter);
808 }
809 findParticleStableDescendants(theDaughter, allstabledescendants);
810 }
811 }
812
816
817 template <class T> bool isHardScatteringVertex(T pVert) {
818 if (pVert == nullptr) return false;
819 T pV = pVert;
820 int numOfPartIn(0);
821 int pdg(0);
822
823 do {
824 pVert = pV;
825 auto incoming = pVert->particles_in();
826 numOfPartIn = incoming.size();
827 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
828 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
829
830 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
831
832 if (numOfPartIn == 2) {
833 auto incoming = pVert->particles_in();
834 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
835 }
836 return false;
837}
838
842
843 template <class T, class U>
844 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
845 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
846 auto vtx = p->production_vertex();
847 if (!vtx) return false;
848 bool fromHad = false;
849 for ( const auto& parent : particles_in(vtx) ) {
850 if (!parent) continue;
851 // should this really go into parton-level territory?
852 // probably depends where BSM particles are being decayed
853 fromBSM |= isBSM(parent);
854 if (!isPhysical(parent)) return false;
855 fromTau |= isTau(parent);
856 if (isHadron(parent)&&!isBeam(parent)) {
857 if (!hadron) hadron = parent; // assumes linear hadron parentage
858 return true;
859 }
860 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
861 }
862 return fromHad;
863 }
864
867
868 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
869 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
870 decltype(thePart->end_vertex()) pVert(nullptr);
871 if (EndVert != nullptr) {
872 do {
873 bool samePart = false;
874 pVert = nullptr;
875 auto outgoing = EndVert->particles_out();
876 auto incoming = EndVert->particles_in();
877 for (const auto& itrDaug: outgoing) {
878 if (!itrDaug) continue;
879 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
880 // brem on generator level for tau
881 (outgoing.size() == 1 && incoming.size() == 1 &&
883 itrDaug->pdg_id() == thePart->pdg_id()) {
884 samePart = true;
885 pVert = itrDaug->end_vertex();
886 }
887 }
888 if (samePart) EndVert = pVert;
889 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
890 }
891 return EndVert;
892 }
893
895
896 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
897 if (!theVert) return {};
898 decltype(theVert->particles_out()) finalStatePart;
899 auto outgoing = theVert->particles_out();
900 for (const auto& thePart: outgoing) {
901 if (!thePart) continue;
902 finalStatePart.push_back(thePart);
903 if (isStable(thePart)) continue;
904 V pVert = findSimulatedEndVertex(thePart);
905 if (pVert == theVert) break; // to prevent Sherpa loop
906 if (pVert != nullptr) {
907 auto vecPart = findFinalStateParticles<V>(pVert);
908 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
909 }
910 }
911 return finalStatePart;
912 }
913#if !defined(XAOD_ANALYSIS)
914#include "AtlasHepMC/GenEvent.h"
915#ifdef HEPMC3
916inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
917inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
918#else
919inline void GeVToMeV(HepMC::GenEvent* evt) {
920 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
921 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
922 (*p)->momentum().py() * 1000,
923 (*p)->momentum().pz() * 1000,
924 (*p)->momentum().e() * 1000);
925 (*p)->set_momentum(fv);
926 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
927 }
928}
929inline void MeVToGeV(HepMC::GenEvent* evt) {
930 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
931 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
932 (*p)->momentum().py() / 1000,
933 (*p)->momentum().pz() / 1000,
934 (*p)->momentum().e() / 1000);
935 (*p)->set_momentum(fv);
936 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
937 }
938}
939#endif
940#endif
941}
942#endif

◆ newGenEvent()

GenEvent * MC::HepMC::newGenEvent ( const int a,
const int b )
inline

Definition at line 624 of file HepMCHelpers.h.

645{
646inline
647auto particles_in (const HepMC::GenVertex* p) {
648 return std::ranges::subrange (p->particles_in_const_begin(),
649 p->particles_in_const_end());
650}
651}
652#endif
653
654namespace MC
655{
656 template <class VTX>
657 auto particles_in (const VTX* p) { return p->particles_in(); }
658 template <class VTX>
659 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
660
661 namespace Pythia8
662 {
664 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
665
666 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
667
668 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
669 }
670
671#include "AtlasPID.h"
672
674 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
675
677 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
678
680 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
681
683 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
684
686 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
687
689 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
690
692 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
693
695 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
696
698 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
699
701 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
702
706 template <class T> inline bool isStableOrSimDecayed(const T& p) {
707 const auto vertex = p->end_vertex();
708 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
709 }
710
712 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
713
715 template <class T> inline bool isSpecialNonInteracting(const T& p) {
716 const int apid = std::abs(p->pdg_id());
717 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
718 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
719 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
720 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
721 return false;
722 }
723
725
726 template <class T> T findMother(T thePart) {
727 auto partOriVert = thePart->production_vertex();
728 if (!partOriVert) return nullptr;
729
730 long partPDG = thePart->pdg_id();
731 long MotherPDG(0);
732
733 auto MothOriVert = partOriVert;
734 MothOriVert = nullptr;
735 T theMoth(nullptr);
736
737 size_t itr = 0;
738 do {
739 if (itr != 0) partOriVert = MothOriVert;
740 for ( const auto& p : particles_in(partOriVert) ) {
741 theMoth = p;
742 if (!theMoth) continue;
743 MotherPDG = theMoth->pdg_id();
744 MothOriVert = theMoth->production_vertex();
745 if (MotherPDG == partPDG) break;
746 }
747 itr++;
748 if (itr > 100) {
749 break;
750 }
751 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
752 MothOriVert != partOriVert);
753 return theMoth;
754 }
755
757
758 template <class C, class T> T findMatching(C TruthContainer, T p) {
759 T ptrPart = nullptr;
760 if (!p) return ptrPart;
761 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
762 for (T truthParticle : *TruthContainer) {
763 if (HepMC::is_sim_descendant(p,truthParticle)) {
764 ptrPart = truthParticle;
765 break;
766 }
767 }
768 }
769 else {
770 for (T truthParticle : TruthContainer) {
771 if (HepMC::is_sim_descendant(p,truthParticle)) {
772 ptrPart = truthParticle;
773 break;
774 }
775 }
776 }
777 return ptrPart;
778 }
780
781 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
782 auto prodVtx = thePart->production_vertex();
783 if (!prodVtx) return;
784 for (const auto& theMother: prodVtx->particles_in()) {
785 if (!theMother) continue;
786 allancestors.insert(theMother);
787 findParticleAncestors(theMother, allancestors);
788 }
789 }
790
792
793 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
794 auto endVtx = thePart->end_vertex();
795 if (!endVtx) return;
796 for (const auto& theDaughter: endVtx->particles_out()) {
797 if (!theDaughter) continue;
798 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
799 allstabledescendants.insert(theDaughter);
800 }
801 findParticleStableDescendants(theDaughter, allstabledescendants);
802 }
803 }
804
808
809 template <class T> bool isHardScatteringVertex(T pVert) {
810 if (pVert == nullptr) return false;
811 T pV = pVert;
812 int numOfPartIn(0);
813 int pdg(0);
814
815 do {
816 pVert = pV;
817 auto incoming = pVert->particles_in();
818 numOfPartIn = incoming.size();
819 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
820 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
821
822 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
823
824 if (numOfPartIn == 2) {
825 auto incoming = pVert->particles_in();
826 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
827 }
828 return false;
829}
830
834
835 template <class T, class U>
836 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
837 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
838 auto vtx = p->production_vertex();
839 if (!vtx) return false;
840 bool fromHad = false;
841 for ( const auto& parent : particles_in(vtx) ) {
842 if (!parent) continue;
843 // should this really go into parton-level territory?
844 // probably depends where BSM particles are being decayed
845 fromBSM |= isBSM(parent);
846 if (!isPhysical(parent)) return false;
847 fromTau |= isTau(parent);
848 if (isHadron(parent)&&!isBeam(parent)) {
849 if (!hadron) hadron = parent; // assumes linear hadron parentage
850 return true;
851 }
852 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
853 }
854 return fromHad;
855 }
856
859
860 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
861 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
862 decltype(thePart->end_vertex()) pVert(nullptr);
863 if (EndVert != nullptr) {
864 do {
865 bool samePart = false;
866 pVert = nullptr;
867 auto outgoing = EndVert->particles_out();
868 auto incoming = EndVert->particles_in();
869 for (const auto& itrDaug: outgoing) {
870 if (!itrDaug) continue;
871 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
872 // brem on generator level for tau
873 (outgoing.size() == 1 && incoming.size() == 1 &&
875 itrDaug->pdg_id() == thePart->pdg_id()) {
876 samePart = true;
877 pVert = itrDaug->end_vertex();
878 }
879 }
880 if (samePart) EndVert = pVert;
881 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
882 }
883 return EndVert;
884 }
885
887
888 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
889 if (!theVert) return {};
890 decltype(theVert->particles_out()) finalStatePart;
891 auto outgoing = theVert->particles_out();
892 for (const auto& thePart: outgoing) {
893 if (!thePart) continue;
894 finalStatePart.push_back(thePart);
895 if (isStable(thePart)) continue;
896 V pVert = findSimulatedEndVertex(thePart);
897 if (pVert == theVert) break; // to prevent Sherpa loop
898 if (pVert != nullptr) {
899 auto vecPart = findFinalStateParticles<V>(pVert);
900 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
901 }
902 }
903 return finalStatePart;
904 }
905#if !defined(XAOD_ANALYSIS)
906#include "AtlasHepMC/GenEvent.h"
907#ifdef HEPMC3
908inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
909inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
910#else
911inline void GeVToMeV(HepMC::GenEvent* evt) {
912 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
913 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
914 (*p)->momentum().py() * 1000,
915 (*p)->momentum().pz() * 1000,
916 (*p)->momentum().e() * 1000);
917 (*p)->set_momentum(fv);
918 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
919 }
920}
921inline void MeVToGeV(HepMC::GenEvent* evt) {
922 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
923 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
924 (*p)->momentum().py() / 1000,
925 (*p)->momentum().pz() / 1000,
926 (*p)->momentum().e() / 1000);
927 (*p)->set_momentum(fv);
928 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
929 }
930}
931#endif
932#endif
933}
934#endif

◆ newGenVertexPtr()

GenVertexPtr MC::HepMC::newGenVertexPtr ( const HepMC::FourVector & pos = HepMC::FourVector(0.0,0.0,0.0,0.0),
const int i = 0 )
inline

Definition at line 64 of file HepMCHelpers.h.

66 { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}

◆ operator<<()

std::ostream & MC::HepMC::operator<< ( std::ostream & os,
const GenVertex * v )
inline

Definition at line 72 of file HepMCHelpers.h.

◆ set_ll_event_number()

bool MC::HepMC::set_ll_event_number ( HepMC::GenEvent * e,
long long int num )
inline

Definition at line 612 of file HepMCHelpers.h.

633 {
634inline
635auto particles_in (const HepMC::GenVertex* p) {
636 return std::ranges::subrange (p->particles_in_const_begin(),
637 p->particles_in_const_end());
638}
639}
640#endif
641
642namespace MC
643{
644 template <class VTX>
645 auto particles_in (const VTX* p) { return p->particles_in(); }
646 template <class VTX>
647 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
648
649 namespace Pythia8
650 {
652 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
653
654 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
655
656 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
657 }
658
659#include "AtlasPID.h"
660
662 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
663
665 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
666
668 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
669
671 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
672
674 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
675
677 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
678
680 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
681
683 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
684
686 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
687
689 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
690
694 template <class T> inline bool isStableOrSimDecayed(const T& p) {
695 const auto vertex = p->end_vertex();
696 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
697 }
698
700 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
701
703 template <class T> inline bool isSpecialNonInteracting(const T& p) {
704 const int apid = std::abs(p->pdg_id());
705 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
706 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
707 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
708 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
709 return false;
710 }
711
713
714 template <class T> T findMother(T thePart) {
715 auto partOriVert = thePart->production_vertex();
716 if (!partOriVert) return nullptr;
717
718 long partPDG = thePart->pdg_id();
719 long MotherPDG(0);
720
721 auto MothOriVert = partOriVert;
722 MothOriVert = nullptr;
723 T theMoth(nullptr);
724
725 size_t itr = 0;
726 do {
727 if (itr != 0) partOriVert = MothOriVert;
728 for ( const auto& p : particles_in(partOriVert) ) {
729 theMoth = p;
730 if (!theMoth) continue;
731 MotherPDG = theMoth->pdg_id();
732 MothOriVert = theMoth->production_vertex();
733 if (MotherPDG == partPDG) break;
734 }
735 itr++;
736 if (itr > 100) {
737 break;
738 }
739 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
740 MothOriVert != partOriVert);
741 return theMoth;
742 }
743
745
746 template <class C, class T> T findMatching(C TruthContainer, T p) {
747 T ptrPart = nullptr;
748 if (!p) return ptrPart;
749 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
750 for (T truthParticle : *TruthContainer) {
751 if (HepMC::is_sim_descendant(p,truthParticle)) {
752 ptrPart = truthParticle;
753 break;
754 }
755 }
756 }
757 else {
758 for (T truthParticle : TruthContainer) {
759 if (HepMC::is_sim_descendant(p,truthParticle)) {
760 ptrPart = truthParticle;
761 break;
762 }
763 }
764 }
765 return ptrPart;
766 }
768
769 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
770 auto prodVtx = thePart->production_vertex();
771 if (!prodVtx) return;
772 for (const auto& theMother: prodVtx->particles_in()) {
773 if (!theMother) continue;
774 allancestors.insert(theMother);
775 findParticleAncestors(theMother, allancestors);
776 }
777 }
778
780
781 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
782 auto endVtx = thePart->end_vertex();
783 if (!endVtx) return;
784 for (const auto& theDaughter: endVtx->particles_out()) {
785 if (!theDaughter) continue;
786 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
787 allstabledescendants.insert(theDaughter);
788 }
789 findParticleStableDescendants(theDaughter, allstabledescendants);
790 }
791 }
792
796
797 template <class T> bool isHardScatteringVertex(T pVert) {
798 if (pVert == nullptr) return false;
799 T pV = pVert;
800 int numOfPartIn(0);
801 int pdg(0);
802
803 do {
804 pVert = pV;
805 auto incoming = pVert->particles_in();
806 numOfPartIn = incoming.size();
807 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
808 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
809
810 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
811
812 if (numOfPartIn == 2) {
813 auto incoming = pVert->particles_in();
814 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
815 }
816 return false;
817}
818
822
823 template <class T, class U>
824 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
825 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
826 auto vtx = p->production_vertex();
827 if (!vtx) return false;
828 bool fromHad = false;
829 for ( const auto& parent : particles_in(vtx) ) {
830 if (!parent) continue;
831 // should this really go into parton-level territory?
832 // probably depends where BSM particles are being decayed
833 fromBSM |= isBSM(parent);
834 if (!isPhysical(parent)) return false;
835 fromTau |= isTau(parent);
836 if (isHadron(parent)&&!isBeam(parent)) {
837 if (!hadron) hadron = parent; // assumes linear hadron parentage
838 return true;
839 }
840 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
841 }
842 return fromHad;
843 }
844
847
848 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
849 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
850 decltype(thePart->end_vertex()) pVert(nullptr);
851 if (EndVert != nullptr) {
852 do {
853 bool samePart = false;
854 pVert = nullptr;
855 auto outgoing = EndVert->particles_out();
856 auto incoming = EndVert->particles_in();
857 for (const auto& itrDaug: outgoing) {
858 if (!itrDaug) continue;
859 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
860 // brem on generator level for tau
861 (outgoing.size() == 1 && incoming.size() == 1 &&
863 itrDaug->pdg_id() == thePart->pdg_id()) {
864 samePart = true;
865 pVert = itrDaug->end_vertex();
866 }
867 }
868 if (samePart) EndVert = pVert;
869 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
870 }
871 return EndVert;
872 }
873
875
876 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
877 if (!theVert) return {};
878 decltype(theVert->particles_out()) finalStatePart;
879 auto outgoing = theVert->particles_out();
880 for (const auto& thePart: outgoing) {
881 if (!thePart) continue;
882 finalStatePart.push_back(thePart);
883 if (isStable(thePart)) continue;
884 V pVert = findSimulatedEndVertex(thePart);
885 if (pVert == theVert) break; // to prevent Sherpa loop
886 if (pVert != nullptr) {
887 auto vecPart = findFinalStateParticles<V>(pVert);
888 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
889 }
890 }
891 return finalStatePart;
892 }
893#if !defined(XAOD_ANALYSIS)
894#include "AtlasHepMC/GenEvent.h"
895#ifdef HEPMC3
896inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
897inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
898#else
899inline void GeVToMeV(HepMC::GenEvent* evt) {
900 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
901 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
902 (*p)->momentum().py() * 1000,
903 (*p)->momentum().pz() * 1000,
904 (*p)->momentum().e() * 1000);
905 (*p)->set_momentum(fv);
906 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
907 }
908}
909inline void MeVToGeV(HepMC::GenEvent* evt) {
910 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
911 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
912 (*p)->momentum().py() / 1000,
913 (*p)->momentum().pz() / 1000,
914 (*p)->momentum().e() / 1000);
915 (*p)->set_momentum(fv);
916 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
917 }
918}
919#endif
920#endif
921}
922#endif

◆ set_mpi()

void MC::HepMC::set_mpi ( GenEvent * e,
const int i )
inline

Definition at line 644 of file HepMCHelpers.h.

665 {
666inline
667auto particles_in (const HepMC::GenVertex* p) {
668 return std::ranges::subrange (p->particles_in_const_begin(),
669 p->particles_in_const_end());
670}
671}
672#endif
673
674namespace MC
675{
676 template <class VTX>
677 auto particles_in (const VTX* p) { return p->particles_in(); }
678 template <class VTX>
679 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
680
681 namespace Pythia8
682 {
684 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
685
686 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
687
688 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
689 }
690
691#include "AtlasPID.h"
692
694 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
695
697 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
698
700 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
701
703 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
704
706 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
707
709 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
710
712 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
713
715 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
716
718 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
719
721 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
722
726 template <class T> inline bool isStableOrSimDecayed(const T& p) {
727 const auto vertex = p->end_vertex();
728 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
729 }
730
732 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
733
735 template <class T> inline bool isSpecialNonInteracting(const T& p) {
736 const int apid = std::abs(p->pdg_id());
737 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
738 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
739 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
740 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
741 return false;
742 }
743
745
746 template <class T> T findMother(T thePart) {
747 auto partOriVert = thePart->production_vertex();
748 if (!partOriVert) return nullptr;
749
750 long partPDG = thePart->pdg_id();
751 long MotherPDG(0);
752
753 auto MothOriVert = partOriVert;
754 MothOriVert = nullptr;
755 T theMoth(nullptr);
756
757 size_t itr = 0;
758 do {
759 if (itr != 0) partOriVert = MothOriVert;
760 for ( const auto& p : particles_in(partOriVert) ) {
761 theMoth = p;
762 if (!theMoth) continue;
763 MotherPDG = theMoth->pdg_id();
764 MothOriVert = theMoth->production_vertex();
765 if (MotherPDG == partPDG) break;
766 }
767 itr++;
768 if (itr > 100) {
769 break;
770 }
771 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
772 MothOriVert != partOriVert);
773 return theMoth;
774 }
775
777
778 template <class C, class T> T findMatching(C TruthContainer, T p) {
779 T ptrPart = nullptr;
780 if (!p) return ptrPart;
781 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
782 for (T truthParticle : *TruthContainer) {
783 if (HepMC::is_sim_descendant(p,truthParticle)) {
784 ptrPart = truthParticle;
785 break;
786 }
787 }
788 }
789 else {
790 for (T truthParticle : TruthContainer) {
791 if (HepMC::is_sim_descendant(p,truthParticle)) {
792 ptrPart = truthParticle;
793 break;
794 }
795 }
796 }
797 return ptrPart;
798 }
800
801 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
802 auto prodVtx = thePart->production_vertex();
803 if (!prodVtx) return;
804 for (const auto& theMother: prodVtx->particles_in()) {
805 if (!theMother) continue;
806 allancestors.insert(theMother);
807 findParticleAncestors(theMother, allancestors);
808 }
809 }
810
812
813 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
814 auto endVtx = thePart->end_vertex();
815 if (!endVtx) return;
816 for (const auto& theDaughter: endVtx->particles_out()) {
817 if (!theDaughter) continue;
818 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
819 allstabledescendants.insert(theDaughter);
820 }
821 findParticleStableDescendants(theDaughter, allstabledescendants);
822 }
823 }
824
828
829 template <class T> bool isHardScatteringVertex(T pVert) {
830 if (pVert == nullptr) return false;
831 T pV = pVert;
832 int numOfPartIn(0);
833 int pdg(0);
834
835 do {
836 pVert = pV;
837 auto incoming = pVert->particles_in();
838 numOfPartIn = incoming.size();
839 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
840 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
841
842 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
843
844 if (numOfPartIn == 2) {
845 auto incoming = pVert->particles_in();
846 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
847 }
848 return false;
849}
850
854
855 template <class T, class U>
856 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
857 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
858 auto vtx = p->production_vertex();
859 if (!vtx) return false;
860 bool fromHad = false;
861 for ( const auto& parent : particles_in(vtx) ) {
862 if (!parent) continue;
863 // should this really go into parton-level territory?
864 // probably depends where BSM particles are being decayed
865 fromBSM |= isBSM(parent);
866 if (!isPhysical(parent)) return false;
867 fromTau |= isTau(parent);
868 if (isHadron(parent)&&!isBeam(parent)) {
869 if (!hadron) hadron = parent; // assumes linear hadron parentage
870 return true;
871 }
872 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
873 }
874 return fromHad;
875 }
876
879
880 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
881 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
882 decltype(thePart->end_vertex()) pVert(nullptr);
883 if (EndVert != nullptr) {
884 do {
885 bool samePart = false;
886 pVert = nullptr;
887 auto outgoing = EndVert->particles_out();
888 auto incoming = EndVert->particles_in();
889 for (const auto& itrDaug: outgoing) {
890 if (!itrDaug) continue;
891 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
892 // brem on generator level for tau
893 (outgoing.size() == 1 && incoming.size() == 1 &&
895 itrDaug->pdg_id() == thePart->pdg_id()) {
896 samePart = true;
897 pVert = itrDaug->end_vertex();
898 }
899 }
900 if (samePart) EndVert = pVert;
901 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
902 }
903 return EndVert;
904 }
905
907
908 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
909 if (!theVert) return {};
910 decltype(theVert->particles_out()) finalStatePart;
911 auto outgoing = theVert->particles_out();
912 for (const auto& thePart: outgoing) {
913 if (!thePart) continue;
914 finalStatePart.push_back(thePart);
915 if (isStable(thePart)) continue;
916 V pVert = findSimulatedEndVertex(thePart);
917 if (pVert == theVert) break; // to prevent Sherpa loop
918 if (pVert != nullptr) {
919 auto vecPart = findFinalStateParticles<V>(pVert);
920 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
921 }
922 }
923 return finalStatePart;
924 }
925#if !defined(XAOD_ANALYSIS)
926#include "AtlasHepMC/GenEvent.h"
927#ifdef HEPMC3
928inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
929inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
930#else
931inline void GeVToMeV(HepMC::GenEvent* evt) {
932 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
933 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
934 (*p)->momentum().py() * 1000,
935 (*p)->momentum().pz() * 1000,
936 (*p)->momentum().e() * 1000);
937 (*p)->set_momentum(fv);
938 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
939 }
940}
941inline void MeVToGeV(HepMC::GenEvent* evt) {
942 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
943 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
944 (*p)->momentum().py() / 1000,
945 (*p)->momentum().pz() / 1000,
946 (*p)->momentum().e() / 1000);
947 (*p)->set_momentum(fv);
948 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
949 }
950}
951#endif
952#endif
953}
954#endif

◆ set_random_states()

template<class T>
void MC::HepMC::set_random_states ( GenEvent * e,
std::vector< T > a )

Definition at line 647 of file HepMCHelpers.h.

668 {
669inline
670auto particles_in (const HepMC::GenVertex* p) {
671 return std::ranges::subrange (p->particles_in_const_begin(),
672 p->particles_in_const_end());
673}
674}
675#endif
676
677namespace MC
678{
679 template <class VTX>
680 auto particles_in (const VTX* p) { return p->particles_in(); }
681 template <class VTX>
682 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
683
684 namespace Pythia8
685 {
687 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
688
689 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
690
691 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
692 }
693
694#include "AtlasPID.h"
695
697 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
698
700 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
701
703 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
704
706 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
707
709 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
710
712 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
713
715 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
716
718 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
719
721 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
722
724 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
725
729 template <class T> inline bool isStableOrSimDecayed(const T& p) {
730 const auto vertex = p->end_vertex();
731 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
732 }
733
735 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
736
738 template <class T> inline bool isSpecialNonInteracting(const T& p) {
739 const int apid = std::abs(p->pdg_id());
740 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
741 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
742 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
743 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
744 return false;
745 }
746
748
749 template <class T> T findMother(T thePart) {
750 auto partOriVert = thePart->production_vertex();
751 if (!partOriVert) return nullptr;
752
753 long partPDG = thePart->pdg_id();
754 long MotherPDG(0);
755
756 auto MothOriVert = partOriVert;
757 MothOriVert = nullptr;
758 T theMoth(nullptr);
759
760 size_t itr = 0;
761 do {
762 if (itr != 0) partOriVert = MothOriVert;
763 for ( const auto& p : particles_in(partOriVert) ) {
764 theMoth = p;
765 if (!theMoth) continue;
766 MotherPDG = theMoth->pdg_id();
767 MothOriVert = theMoth->production_vertex();
768 if (MotherPDG == partPDG) break;
769 }
770 itr++;
771 if (itr > 100) {
772 break;
773 }
774 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
775 MothOriVert != partOriVert);
776 return theMoth;
777 }
778
780
781 template <class C, class T> T findMatching(C TruthContainer, T p) {
782 T ptrPart = nullptr;
783 if (!p) return ptrPart;
784 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
785 for (T truthParticle : *TruthContainer) {
786 if (HepMC::is_sim_descendant(p,truthParticle)) {
787 ptrPart = truthParticle;
788 break;
789 }
790 }
791 }
792 else {
793 for (T truthParticle : TruthContainer) {
794 if (HepMC::is_sim_descendant(p,truthParticle)) {
795 ptrPart = truthParticle;
796 break;
797 }
798 }
799 }
800 return ptrPart;
801 }
803
804 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
805 auto prodVtx = thePart->production_vertex();
806 if (!prodVtx) return;
807 for (const auto& theMother: prodVtx->particles_in()) {
808 if (!theMother) continue;
809 allancestors.insert(theMother);
810 findParticleAncestors(theMother, allancestors);
811 }
812 }
813
815
816 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
817 auto endVtx = thePart->end_vertex();
818 if (!endVtx) return;
819 for (const auto& theDaughter: endVtx->particles_out()) {
820 if (!theDaughter) continue;
821 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
822 allstabledescendants.insert(theDaughter);
823 }
824 findParticleStableDescendants(theDaughter, allstabledescendants);
825 }
826 }
827
831
832 template <class T> bool isHardScatteringVertex(T pVert) {
833 if (pVert == nullptr) return false;
834 T pV = pVert;
835 int numOfPartIn(0);
836 int pdg(0);
837
838 do {
839 pVert = pV;
840 auto incoming = pVert->particles_in();
841 numOfPartIn = incoming.size();
842 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
843 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
844
845 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
846
847 if (numOfPartIn == 2) {
848 auto incoming = pVert->particles_in();
849 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
850 }
851 return false;
852}
853
857
858 template <class T, class U>
859 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
860 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
861 auto vtx = p->production_vertex();
862 if (!vtx) return false;
863 bool fromHad = false;
864 for ( const auto& parent : particles_in(vtx) ) {
865 if (!parent) continue;
866 // should this really go into parton-level territory?
867 // probably depends where BSM particles are being decayed
868 fromBSM |= isBSM(parent);
869 if (!isPhysical(parent)) return false;
870 fromTau |= isTau(parent);
871 if (isHadron(parent)&&!isBeam(parent)) {
872 if (!hadron) hadron = parent; // assumes linear hadron parentage
873 return true;
874 }
875 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
876 }
877 return fromHad;
878 }
879
882
883 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
884 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
885 decltype(thePart->end_vertex()) pVert(nullptr);
886 if (EndVert != nullptr) {
887 do {
888 bool samePart = false;
889 pVert = nullptr;
890 auto outgoing = EndVert->particles_out();
891 auto incoming = EndVert->particles_in();
892 for (const auto& itrDaug: outgoing) {
893 if (!itrDaug) continue;
894 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
895 // brem on generator level for tau
896 (outgoing.size() == 1 && incoming.size() == 1 &&
898 itrDaug->pdg_id() == thePart->pdg_id()) {
899 samePart = true;
900 pVert = itrDaug->end_vertex();
901 }
902 }
903 if (samePart) EndVert = pVert;
904 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
905 }
906 return EndVert;
907 }
908
910
911 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
912 if (!theVert) return {};
913 decltype(theVert->particles_out()) finalStatePart;
914 auto outgoing = theVert->particles_out();
915 for (const auto& thePart: outgoing) {
916 if (!thePart) continue;
917 finalStatePart.push_back(thePart);
918 if (isStable(thePart)) continue;
919 V pVert = findSimulatedEndVertex(thePart);
920 if (pVert == theVert) break; // to prevent Sherpa loop
921 if (pVert != nullptr) {
922 auto vecPart = findFinalStateParticles<V>(pVert);
923 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
924 }
925 }
926 return finalStatePart;
927 }
928#if !defined(XAOD_ANALYSIS)
929#include "AtlasHepMC/GenEvent.h"
930#ifdef HEPMC3
931inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
932inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
933#else
934inline void GeVToMeV(HepMC::GenEvent* evt) {
935 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
936 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
937 (*p)->momentum().py() * 1000,
938 (*p)->momentum().pz() * 1000,
939 (*p)->momentum().e() * 1000);
940 (*p)->set_momentum(fv);
941 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
942 }
943}
944inline void MeVToGeV(HepMC::GenEvent* evt) {
945 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
946 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
947 (*p)->momentum().py() / 1000,
948 (*p)->momentum().pz() / 1000,
949 (*p)->momentum().e() / 1000);
950 (*p)->set_momentum(fv);
951 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
952 }
953}
954#endif
955#endif
956}
957#endif

◆ set_signal_process_id()

void MC::HepMC::set_signal_process_id ( GenEvent * e,
const int i )
inline

Definition at line 641 of file HepMCHelpers.h.

662 {
663inline
664auto particles_in (const HepMC::GenVertex* p) {
665 return std::ranges::subrange (p->particles_in_const_begin(),
666 p->particles_in_const_end());
667}
668}
669#endif
670
671namespace MC
672{
673 template <class VTX>
674 auto particles_in (const VTX* p) { return p->particles_in(); }
675 template <class VTX>
676 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
677
678 namespace Pythia8
679 {
681 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
682
683 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
684
685 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
686 }
687
688#include "AtlasPID.h"
689
691 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
692
694 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
695
697 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
698
700 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
701
703 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
704
706 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
707
709 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
710
712 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
713
715 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
716
718 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
719
723 template <class T> inline bool isStableOrSimDecayed(const T& p) {
724 const auto vertex = p->end_vertex();
725 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
726 }
727
729 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
730
732 template <class T> inline bool isSpecialNonInteracting(const T& p) {
733 const int apid = std::abs(p->pdg_id());
734 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
735 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
736 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
737 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
738 return false;
739 }
740
742
743 template <class T> T findMother(T thePart) {
744 auto partOriVert = thePart->production_vertex();
745 if (!partOriVert) return nullptr;
746
747 long partPDG = thePart->pdg_id();
748 long MotherPDG(0);
749
750 auto MothOriVert = partOriVert;
751 MothOriVert = nullptr;
752 T theMoth(nullptr);
753
754 size_t itr = 0;
755 do {
756 if (itr != 0) partOriVert = MothOriVert;
757 for ( const auto& p : particles_in(partOriVert) ) {
758 theMoth = p;
759 if (!theMoth) continue;
760 MotherPDG = theMoth->pdg_id();
761 MothOriVert = theMoth->production_vertex();
762 if (MotherPDG == partPDG) break;
763 }
764 itr++;
765 if (itr > 100) {
766 break;
767 }
768 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
769 MothOriVert != partOriVert);
770 return theMoth;
771 }
772
774
775 template <class C, class T> T findMatching(C TruthContainer, T p) {
776 T ptrPart = nullptr;
777 if (!p) return ptrPart;
778 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
779 for (T truthParticle : *TruthContainer) {
780 if (HepMC::is_sim_descendant(p,truthParticle)) {
781 ptrPart = truthParticle;
782 break;
783 }
784 }
785 }
786 else {
787 for (T truthParticle : TruthContainer) {
788 if (HepMC::is_sim_descendant(p,truthParticle)) {
789 ptrPart = truthParticle;
790 break;
791 }
792 }
793 }
794 return ptrPart;
795 }
797
798 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
799 auto prodVtx = thePart->production_vertex();
800 if (!prodVtx) return;
801 for (const auto& theMother: prodVtx->particles_in()) {
802 if (!theMother) continue;
803 allancestors.insert(theMother);
804 findParticleAncestors(theMother, allancestors);
805 }
806 }
807
809
810 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
811 auto endVtx = thePart->end_vertex();
812 if (!endVtx) return;
813 for (const auto& theDaughter: endVtx->particles_out()) {
814 if (!theDaughter) continue;
815 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
816 allstabledescendants.insert(theDaughter);
817 }
818 findParticleStableDescendants(theDaughter, allstabledescendants);
819 }
820 }
821
825
826 template <class T> bool isHardScatteringVertex(T pVert) {
827 if (pVert == nullptr) return false;
828 T pV = pVert;
829 int numOfPartIn(0);
830 int pdg(0);
831
832 do {
833 pVert = pV;
834 auto incoming = pVert->particles_in();
835 numOfPartIn = incoming.size();
836 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
837 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
838
839 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
840
841 if (numOfPartIn == 2) {
842 auto incoming = pVert->particles_in();
843 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
844 }
845 return false;
846}
847
851
852 template <class T, class U>
853 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
854 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
855 auto vtx = p->production_vertex();
856 if (!vtx) return false;
857 bool fromHad = false;
858 for ( const auto& parent : particles_in(vtx) ) {
859 if (!parent) continue;
860 // should this really go into parton-level territory?
861 // probably depends where BSM particles are being decayed
862 fromBSM |= isBSM(parent);
863 if (!isPhysical(parent)) return false;
864 fromTau |= isTau(parent);
865 if (isHadron(parent)&&!isBeam(parent)) {
866 if (!hadron) hadron = parent; // assumes linear hadron parentage
867 return true;
868 }
869 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
870 }
871 return fromHad;
872 }
873
876
877 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
878 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
879 decltype(thePart->end_vertex()) pVert(nullptr);
880 if (EndVert != nullptr) {
881 do {
882 bool samePart = false;
883 pVert = nullptr;
884 auto outgoing = EndVert->particles_out();
885 auto incoming = EndVert->particles_in();
886 for (const auto& itrDaug: outgoing) {
887 if (!itrDaug) continue;
888 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
889 // brem on generator level for tau
890 (outgoing.size() == 1 && incoming.size() == 1 &&
892 itrDaug->pdg_id() == thePart->pdg_id()) {
893 samePart = true;
894 pVert = itrDaug->end_vertex();
895 }
896 }
897 if (samePart) EndVert = pVert;
898 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
899 }
900 return EndVert;
901 }
902
904
905 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
906 if (!theVert) return {};
907 decltype(theVert->particles_out()) finalStatePart;
908 auto outgoing = theVert->particles_out();
909 for (const auto& thePart: outgoing) {
910 if (!thePart) continue;
911 finalStatePart.push_back(thePart);
912 if (isStable(thePart)) continue;
913 V pVert = findSimulatedEndVertex(thePart);
914 if (pVert == theVert) break; // to prevent Sherpa loop
915 if (pVert != nullptr) {
916 auto vecPart = findFinalStateParticles<V>(pVert);
917 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
918 }
919 }
920 return finalStatePart;
921 }
922#if !defined(XAOD_ANALYSIS)
923#include "AtlasHepMC/GenEvent.h"
924#ifdef HEPMC3
925inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
926inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
927#else
928inline void GeVToMeV(HepMC::GenEvent* evt) {
929 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
930 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
931 (*p)->momentum().py() * 1000,
932 (*p)->momentum().pz() * 1000,
933 (*p)->momentum().e() * 1000);
934 (*p)->set_momentum(fv);
935 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
936 }
937}
938inline void MeVToGeV(HepMC::GenEvent* evt) {
939 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
940 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
941 (*p)->momentum().py() / 1000,
942 (*p)->momentum().pz() / 1000,
943 (*p)->momentum().e() / 1000);
944 (*p)->set_momentum(fv);
945 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
946 }
947}
948#endif
949#endif
950}
951#endif

◆ set_signal_process_vertex()

template<class T>
void MC::HepMC::set_signal_process_vertex ( GenEvent * e,
T v )

Definition at line 650 of file HepMCHelpers.h.

671 {
672inline
673auto particles_in (const HepMC::GenVertex* p) {
674 return std::ranges::subrange (p->particles_in_const_begin(),
675 p->particles_in_const_end());
676}
677}
678#endif
679
680namespace MC
681{
682 template <class VTX>
683 auto particles_in (const VTX* p) { return p->particles_in(); }
684 template <class VTX>
685 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
686
687 namespace Pythia8
688 {
690 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
691
692 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
693
694 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
695 }
696
697#include "AtlasPID.h"
698
700 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
701
703 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
704
706 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
707
709 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
710
712 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
713
715 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
716
718 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
719
721 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
722
724 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
725
727 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
728
732 template <class T> inline bool isStableOrSimDecayed(const T& p) {
733 const auto vertex = p->end_vertex();
734 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
735 }
736
738 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
739
741 template <class T> inline bool isSpecialNonInteracting(const T& p) {
742 const int apid = std::abs(p->pdg_id());
743 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
744 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
745 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
746 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
747 return false;
748 }
749
751
752 template <class T> T findMother(T thePart) {
753 auto partOriVert = thePart->production_vertex();
754 if (!partOriVert) return nullptr;
755
756 long partPDG = thePart->pdg_id();
757 long MotherPDG(0);
758
759 auto MothOriVert = partOriVert;
760 MothOriVert = nullptr;
761 T theMoth(nullptr);
762
763 size_t itr = 0;
764 do {
765 if (itr != 0) partOriVert = MothOriVert;
766 for ( const auto& p : particles_in(partOriVert) ) {
767 theMoth = p;
768 if (!theMoth) continue;
769 MotherPDG = theMoth->pdg_id();
770 MothOriVert = theMoth->production_vertex();
771 if (MotherPDG == partPDG) break;
772 }
773 itr++;
774 if (itr > 100) {
775 break;
776 }
777 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
778 MothOriVert != partOriVert);
779 return theMoth;
780 }
781
783
784 template <class C, class T> T findMatching(C TruthContainer, T p) {
785 T ptrPart = nullptr;
786 if (!p) return ptrPart;
787 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
788 for (T truthParticle : *TruthContainer) {
789 if (HepMC::is_sim_descendant(p,truthParticle)) {
790 ptrPart = truthParticle;
791 break;
792 }
793 }
794 }
795 else {
796 for (T truthParticle : TruthContainer) {
797 if (HepMC::is_sim_descendant(p,truthParticle)) {
798 ptrPart = truthParticle;
799 break;
800 }
801 }
802 }
803 return ptrPart;
804 }
806
807 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
808 auto prodVtx = thePart->production_vertex();
809 if (!prodVtx) return;
810 for (const auto& theMother: prodVtx->particles_in()) {
811 if (!theMother) continue;
812 allancestors.insert(theMother);
813 findParticleAncestors(theMother, allancestors);
814 }
815 }
816
818
819 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
820 auto endVtx = thePart->end_vertex();
821 if (!endVtx) return;
822 for (const auto& theDaughter: endVtx->particles_out()) {
823 if (!theDaughter) continue;
824 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
825 allstabledescendants.insert(theDaughter);
826 }
827 findParticleStableDescendants(theDaughter, allstabledescendants);
828 }
829 }
830
834
835 template <class T> bool isHardScatteringVertex(T pVert) {
836 if (pVert == nullptr) return false;
837 T pV = pVert;
838 int numOfPartIn(0);
839 int pdg(0);
840
841 do {
842 pVert = pV;
843 auto incoming = pVert->particles_in();
844 numOfPartIn = incoming.size();
845 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
846 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
847
848 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
849
850 if (numOfPartIn == 2) {
851 auto incoming = pVert->particles_in();
852 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
853 }
854 return false;
855}
856
860
861 template <class T, class U>
862 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
863 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
864 auto vtx = p->production_vertex();
865 if (!vtx) return false;
866 bool fromHad = false;
867 for ( const auto& parent : particles_in(vtx) ) {
868 if (!parent) continue;
869 // should this really go into parton-level territory?
870 // probably depends where BSM particles are being decayed
871 fromBSM |= isBSM(parent);
872 if (!isPhysical(parent)) return false;
873 fromTau |= isTau(parent);
874 if (isHadron(parent)&&!isBeam(parent)) {
875 if (!hadron) hadron = parent; // assumes linear hadron parentage
876 return true;
877 }
878 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
879 }
880 return fromHad;
881 }
882
885
886 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
887 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
888 decltype(thePart->end_vertex()) pVert(nullptr);
889 if (EndVert != nullptr) {
890 do {
891 bool samePart = false;
892 pVert = nullptr;
893 auto outgoing = EndVert->particles_out();
894 auto incoming = EndVert->particles_in();
895 for (const auto& itrDaug: outgoing) {
896 if (!itrDaug) continue;
897 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
898 // brem on generator level for tau
899 (outgoing.size() == 1 && incoming.size() == 1 &&
901 itrDaug->pdg_id() == thePart->pdg_id()) {
902 samePart = true;
903 pVert = itrDaug->end_vertex();
904 }
905 }
906 if (samePart) EndVert = pVert;
907 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
908 }
909 return EndVert;
910 }
911
913
914 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
915 if (!theVert) return {};
916 decltype(theVert->particles_out()) finalStatePart;
917 auto outgoing = theVert->particles_out();
918 for (const auto& thePart: outgoing) {
919 if (!thePart) continue;
920 finalStatePart.push_back(thePart);
921 if (isStable(thePart)) continue;
922 V pVert = findSimulatedEndVertex(thePart);
923 if (pVert == theVert) break; // to prevent Sherpa loop
924 if (pVert != nullptr) {
925 auto vecPart = findFinalStateParticles<V>(pVert);
926 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
927 }
928 }
929 return finalStatePart;
930 }
931#if !defined(XAOD_ANALYSIS)
932#include "AtlasHepMC/GenEvent.h"
933#ifdef HEPMC3
934inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
935inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
936#else
937inline void GeVToMeV(HepMC::GenEvent* evt) {
938 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
939 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
940 (*p)->momentum().py() * 1000,
941 (*p)->momentum().pz() * 1000,
942 (*p)->momentum().e() * 1000);
943 (*p)->set_momentum(fv);
944 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
945 }
946}
947inline void MeVToGeV(HepMC::GenEvent* evt) {
948 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
949 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
950 (*p)->momentum().py() / 1000,
951 (*p)->momentum().pz() / 1000,
952 (*p)->momentum().e() / 1000);
953 (*p)->set_momentum(fv);
954 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
955 }
956}
957#endif
958#endif
959}
960#endif

◆ signal_process_id() [1/2]

int MC::HepMC::signal_process_id ( const GenEvent & e)
inline

Definition at line 635 of file HepMCHelpers.h.

656 {
657inline
658auto particles_in (const HepMC::GenVertex* p) {
659 return std::ranges::subrange (p->particles_in_const_begin(),
660 p->particles_in_const_end());
661}
662}
663#endif
664
665namespace MC
666{
667 template <class VTX>
668 auto particles_in (const VTX* p) { return p->particles_in(); }
669 template <class VTX>
670 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
671
672 namespace Pythia8
673 {
675 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
676
677 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
678
679 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
680 }
681
682#include "AtlasPID.h"
683
685 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
686
688 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
689
691 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
692
694 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
695
697 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
698
700 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
701
703 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
704
706 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
707
709 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
710
712 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
713
717 template <class T> inline bool isStableOrSimDecayed(const T& p) {
718 const auto vertex = p->end_vertex();
719 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
720 }
721
723 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
724
726 template <class T> inline bool isSpecialNonInteracting(const T& p) {
727 const int apid = std::abs(p->pdg_id());
728 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
729 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
730 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
731 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
732 return false;
733 }
734
736
737 template <class T> T findMother(T thePart) {
738 auto partOriVert = thePart->production_vertex();
739 if (!partOriVert) return nullptr;
740
741 long partPDG = thePart->pdg_id();
742 long MotherPDG(0);
743
744 auto MothOriVert = partOriVert;
745 MothOriVert = nullptr;
746 T theMoth(nullptr);
747
748 size_t itr = 0;
749 do {
750 if (itr != 0) partOriVert = MothOriVert;
751 for ( const auto& p : particles_in(partOriVert) ) {
752 theMoth = p;
753 if (!theMoth) continue;
754 MotherPDG = theMoth->pdg_id();
755 MothOriVert = theMoth->production_vertex();
756 if (MotherPDG == partPDG) break;
757 }
758 itr++;
759 if (itr > 100) {
760 break;
761 }
762 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
763 MothOriVert != partOriVert);
764 return theMoth;
765 }
766
768
769 template <class C, class T> T findMatching(C TruthContainer, T p) {
770 T ptrPart = nullptr;
771 if (!p) return ptrPart;
772 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
773 for (T truthParticle : *TruthContainer) {
774 if (HepMC::is_sim_descendant(p,truthParticle)) {
775 ptrPart = truthParticle;
776 break;
777 }
778 }
779 }
780 else {
781 for (T truthParticle : TruthContainer) {
782 if (HepMC::is_sim_descendant(p,truthParticle)) {
783 ptrPart = truthParticle;
784 break;
785 }
786 }
787 }
788 return ptrPart;
789 }
791
792 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
793 auto prodVtx = thePart->production_vertex();
794 if (!prodVtx) return;
795 for (const auto& theMother: prodVtx->particles_in()) {
796 if (!theMother) continue;
797 allancestors.insert(theMother);
798 findParticleAncestors(theMother, allancestors);
799 }
800 }
801
803
804 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
805 auto endVtx = thePart->end_vertex();
806 if (!endVtx) return;
807 for (const auto& theDaughter: endVtx->particles_out()) {
808 if (!theDaughter) continue;
809 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
810 allstabledescendants.insert(theDaughter);
811 }
812 findParticleStableDescendants(theDaughter, allstabledescendants);
813 }
814 }
815
819
820 template <class T> bool isHardScatteringVertex(T pVert) {
821 if (pVert == nullptr) return false;
822 T pV = pVert;
823 int numOfPartIn(0);
824 int pdg(0);
825
826 do {
827 pVert = pV;
828 auto incoming = pVert->particles_in();
829 numOfPartIn = incoming.size();
830 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
831 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
832
833 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
834
835 if (numOfPartIn == 2) {
836 auto incoming = pVert->particles_in();
837 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
838 }
839 return false;
840}
841
845
846 template <class T, class U>
847 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
848 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
849 auto vtx = p->production_vertex();
850 if (!vtx) return false;
851 bool fromHad = false;
852 for ( const auto& parent : particles_in(vtx) ) {
853 if (!parent) continue;
854 // should this really go into parton-level territory?
855 // probably depends where BSM particles are being decayed
856 fromBSM |= isBSM(parent);
857 if (!isPhysical(parent)) return false;
858 fromTau |= isTau(parent);
859 if (isHadron(parent)&&!isBeam(parent)) {
860 if (!hadron) hadron = parent; // assumes linear hadron parentage
861 return true;
862 }
863 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
864 }
865 return fromHad;
866 }
867
870
871 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
872 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
873 decltype(thePart->end_vertex()) pVert(nullptr);
874 if (EndVert != nullptr) {
875 do {
876 bool samePart = false;
877 pVert = nullptr;
878 auto outgoing = EndVert->particles_out();
879 auto incoming = EndVert->particles_in();
880 for (const auto& itrDaug: outgoing) {
881 if (!itrDaug) continue;
882 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
883 // brem on generator level for tau
884 (outgoing.size() == 1 && incoming.size() == 1 &&
886 itrDaug->pdg_id() == thePart->pdg_id()) {
887 samePart = true;
888 pVert = itrDaug->end_vertex();
889 }
890 }
891 if (samePart) EndVert = pVert;
892 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
893 }
894 return EndVert;
895 }
896
898
899 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
900 if (!theVert) return {};
901 decltype(theVert->particles_out()) finalStatePart;
902 auto outgoing = theVert->particles_out();
903 for (const auto& thePart: outgoing) {
904 if (!thePart) continue;
905 finalStatePart.push_back(thePart);
906 if (isStable(thePart)) continue;
907 V pVert = findSimulatedEndVertex(thePart);
908 if (pVert == theVert) break; // to prevent Sherpa loop
909 if (pVert != nullptr) {
910 auto vecPart = findFinalStateParticles<V>(pVert);
911 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
912 }
913 }
914 return finalStatePart;
915 }
916#if !defined(XAOD_ANALYSIS)
917#include "AtlasHepMC/GenEvent.h"
918#ifdef HEPMC3
919inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
920inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
921#else
922inline void GeVToMeV(HepMC::GenEvent* evt) {
923 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
924 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
925 (*p)->momentum().py() * 1000,
926 (*p)->momentum().pz() * 1000,
927 (*p)->momentum().e() * 1000);
928 (*p)->set_momentum(fv);
929 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
930 }
931}
932inline void MeVToGeV(HepMC::GenEvent* evt) {
933 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
934 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
935 (*p)->momentum().py() / 1000,
936 (*p)->momentum().pz() / 1000,
937 (*p)->momentum().e() / 1000);
938 (*p)->set_momentum(fv);
939 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
940 }
941}
942#endif
943#endif
944}
945#endif

◆ signal_process_id() [2/2]

int MC::HepMC::signal_process_id ( const GenEvent * e)
inline

Definition at line 638 of file HepMCHelpers.h.

659 {
660inline
661auto particles_in (const HepMC::GenVertex* p) {
662 return std::ranges::subrange (p->particles_in_const_begin(),
663 p->particles_in_const_end());
664}
665}
666#endif
667
668namespace MC
669{
670 template <class VTX>
671 auto particles_in (const VTX* p) { return p->particles_in(); }
672 template <class VTX>
673 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
674
675 namespace Pythia8
676 {
678 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
679
680 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
681
682 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
683 }
684
685#include "AtlasPID.h"
686
688 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
689
691 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
692
694 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
695
697 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
698
700 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
701
703 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
704
706 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
707
709 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
710
712 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
713
715 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
716
720 template <class T> inline bool isStableOrSimDecayed(const T& p) {
721 const auto vertex = p->end_vertex();
722 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
723 }
724
726 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
727
729 template <class T> inline bool isSpecialNonInteracting(const T& p) {
730 const int apid = std::abs(p->pdg_id());
731 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
732 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
733 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
734 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
735 return false;
736 }
737
739
740 template <class T> T findMother(T thePart) {
741 auto partOriVert = thePart->production_vertex();
742 if (!partOriVert) return nullptr;
743
744 long partPDG = thePart->pdg_id();
745 long MotherPDG(0);
746
747 auto MothOriVert = partOriVert;
748 MothOriVert = nullptr;
749 T theMoth(nullptr);
750
751 size_t itr = 0;
752 do {
753 if (itr != 0) partOriVert = MothOriVert;
754 for ( const auto& p : particles_in(partOriVert) ) {
755 theMoth = p;
756 if (!theMoth) continue;
757 MotherPDG = theMoth->pdg_id();
758 MothOriVert = theMoth->production_vertex();
759 if (MotherPDG == partPDG) break;
760 }
761 itr++;
762 if (itr > 100) {
763 break;
764 }
765 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
766 MothOriVert != partOriVert);
767 return theMoth;
768 }
769
771
772 template <class C, class T> T findMatching(C TruthContainer, T p) {
773 T ptrPart = nullptr;
774 if (!p) return ptrPart;
775 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
776 for (T truthParticle : *TruthContainer) {
777 if (HepMC::is_sim_descendant(p,truthParticle)) {
778 ptrPart = truthParticle;
779 break;
780 }
781 }
782 }
783 else {
784 for (T truthParticle : TruthContainer) {
785 if (HepMC::is_sim_descendant(p,truthParticle)) {
786 ptrPart = truthParticle;
787 break;
788 }
789 }
790 }
791 return ptrPart;
792 }
794
795 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
796 auto prodVtx = thePart->production_vertex();
797 if (!prodVtx) return;
798 for (const auto& theMother: prodVtx->particles_in()) {
799 if (!theMother) continue;
800 allancestors.insert(theMother);
801 findParticleAncestors(theMother, allancestors);
802 }
803 }
804
806
807 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
808 auto endVtx = thePart->end_vertex();
809 if (!endVtx) return;
810 for (const auto& theDaughter: endVtx->particles_out()) {
811 if (!theDaughter) continue;
812 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
813 allstabledescendants.insert(theDaughter);
814 }
815 findParticleStableDescendants(theDaughter, allstabledescendants);
816 }
817 }
818
822
823 template <class T> bool isHardScatteringVertex(T pVert) {
824 if (pVert == nullptr) return false;
825 T pV = pVert;
826 int numOfPartIn(0);
827 int pdg(0);
828
829 do {
830 pVert = pV;
831 auto incoming = pVert->particles_in();
832 numOfPartIn = incoming.size();
833 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
834 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
835
836 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
837
838 if (numOfPartIn == 2) {
839 auto incoming = pVert->particles_in();
840 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
841 }
842 return false;
843}
844
848
849 template <class T, class U>
850 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
851 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
852 auto vtx = p->production_vertex();
853 if (!vtx) return false;
854 bool fromHad = false;
855 for ( const auto& parent : particles_in(vtx) ) {
856 if (!parent) continue;
857 // should this really go into parton-level territory?
858 // probably depends where BSM particles are being decayed
859 fromBSM |= isBSM(parent);
860 if (!isPhysical(parent)) return false;
861 fromTau |= isTau(parent);
862 if (isHadron(parent)&&!isBeam(parent)) {
863 if (!hadron) hadron = parent; // assumes linear hadron parentage
864 return true;
865 }
866 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
867 }
868 return fromHad;
869 }
870
873
874 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
875 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
876 decltype(thePart->end_vertex()) pVert(nullptr);
877 if (EndVert != nullptr) {
878 do {
879 bool samePart = false;
880 pVert = nullptr;
881 auto outgoing = EndVert->particles_out();
882 auto incoming = EndVert->particles_in();
883 for (const auto& itrDaug: outgoing) {
884 if (!itrDaug) continue;
885 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
886 // brem on generator level for tau
887 (outgoing.size() == 1 && incoming.size() == 1 &&
889 itrDaug->pdg_id() == thePart->pdg_id()) {
890 samePart = true;
891 pVert = itrDaug->end_vertex();
892 }
893 }
894 if (samePart) EndVert = pVert;
895 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
896 }
897 return EndVert;
898 }
899
901
902 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
903 if (!theVert) return {};
904 decltype(theVert->particles_out()) finalStatePart;
905 auto outgoing = theVert->particles_out();
906 for (const auto& thePart: outgoing) {
907 if (!thePart) continue;
908 finalStatePart.push_back(thePart);
909 if (isStable(thePart)) continue;
910 V pVert = findSimulatedEndVertex(thePart);
911 if (pVert == theVert) break; // to prevent Sherpa loop
912 if (pVert != nullptr) {
913 auto vecPart = findFinalStateParticles<V>(pVert);
914 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
915 }
916 }
917 return finalStatePart;
918 }
919#if !defined(XAOD_ANALYSIS)
920#include "AtlasHepMC/GenEvent.h"
921#ifdef HEPMC3
922inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
923inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
924#else
925inline void GeVToMeV(HepMC::GenEvent* evt) {
926 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
927 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
928 (*p)->momentum().py() * 1000,
929 (*p)->momentum().pz() * 1000,
930 (*p)->momentum().e() * 1000);
931 (*p)->set_momentum(fv);
932 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
933 }
934}
935inline void MeVToGeV(HepMC::GenEvent* evt) {
936 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
937 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
938 (*p)->momentum().py() / 1000,
939 (*p)->momentum().pz() / 1000,
940 (*p)->momentum().e() / 1000);
941 (*p)->set_momentum(fv);
942 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
943 }
944}
945#endif
946#endif
947}
948#endif

◆ signal_process_vertex()

GenVertex * MC::HepMC::signal_process_vertex ( const GenEvent * e)
inline

Definition at line 625 of file HepMCHelpers.h.

646{
647inline
648auto particles_in (const HepMC::GenVertex* p) {
649 return std::ranges::subrange (p->particles_in_const_begin(),
650 p->particles_in_const_end());
651}
652}
653#endif
654
655namespace MC
656{
657 template <class VTX>
658 auto particles_in (const VTX* p) { return p->particles_in(); }
659 template <class VTX>
660 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
661
662 namespace Pythia8
663 {
665 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
666
667 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
668
669 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
670 }
671
672#include "AtlasPID.h"
673
675 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
676
678 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
679
681 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
682
684 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
685
687 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
688
690 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
691
693 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
694
696 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
697
699 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
700
702 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
703
707 template <class T> inline bool isStableOrSimDecayed(const T& p) {
708 const auto vertex = p->end_vertex();
709 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
710 }
711
713 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
714
716 template <class T> inline bool isSpecialNonInteracting(const T& p) {
717 const int apid = std::abs(p->pdg_id());
718 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
719 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
720 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
721 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
722 return false;
723 }
724
726
727 template <class T> T findMother(T thePart) {
728 auto partOriVert = thePart->production_vertex();
729 if (!partOriVert) return nullptr;
730
731 long partPDG = thePart->pdg_id();
732 long MotherPDG(0);
733
734 auto MothOriVert = partOriVert;
735 MothOriVert = nullptr;
736 T theMoth(nullptr);
737
738 size_t itr = 0;
739 do {
740 if (itr != 0) partOriVert = MothOriVert;
741 for ( const auto& p : particles_in(partOriVert) ) {
742 theMoth = p;
743 if (!theMoth) continue;
744 MotherPDG = theMoth->pdg_id();
745 MothOriVert = theMoth->production_vertex();
746 if (MotherPDG == partPDG) break;
747 }
748 itr++;
749 if (itr > 100) {
750 break;
751 }
752 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
753 MothOriVert != partOriVert);
754 return theMoth;
755 }
756
758
759 template <class C, class T> T findMatching(C TruthContainer, T p) {
760 T ptrPart = nullptr;
761 if (!p) return ptrPart;
762 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
763 for (T truthParticle : *TruthContainer) {
764 if (HepMC::is_sim_descendant(p,truthParticle)) {
765 ptrPart = truthParticle;
766 break;
767 }
768 }
769 }
770 else {
771 for (T truthParticle : TruthContainer) {
772 if (HepMC::is_sim_descendant(p,truthParticle)) {
773 ptrPart = truthParticle;
774 break;
775 }
776 }
777 }
778 return ptrPart;
779 }
781
782 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
783 auto prodVtx = thePart->production_vertex();
784 if (!prodVtx) return;
785 for (const auto& theMother: prodVtx->particles_in()) {
786 if (!theMother) continue;
787 allancestors.insert(theMother);
788 findParticleAncestors(theMother, allancestors);
789 }
790 }
791
793
794 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
795 auto endVtx = thePart->end_vertex();
796 if (!endVtx) return;
797 for (const auto& theDaughter: endVtx->particles_out()) {
798 if (!theDaughter) continue;
799 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
800 allstabledescendants.insert(theDaughter);
801 }
802 findParticleStableDescendants(theDaughter, allstabledescendants);
803 }
804 }
805
809
810 template <class T> bool isHardScatteringVertex(T pVert) {
811 if (pVert == nullptr) return false;
812 T pV = pVert;
813 int numOfPartIn(0);
814 int pdg(0);
815
816 do {
817 pVert = pV;
818 auto incoming = pVert->particles_in();
819 numOfPartIn = incoming.size();
820 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
821 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
822
823 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
824
825 if (numOfPartIn == 2) {
826 auto incoming = pVert->particles_in();
827 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
828 }
829 return false;
830}
831
835
836 template <class T, class U>
837 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
838 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
839 auto vtx = p->production_vertex();
840 if (!vtx) return false;
841 bool fromHad = false;
842 for ( const auto& parent : particles_in(vtx) ) {
843 if (!parent) continue;
844 // should this really go into parton-level territory?
845 // probably depends where BSM particles are being decayed
846 fromBSM |= isBSM(parent);
847 if (!isPhysical(parent)) return false;
848 fromTau |= isTau(parent);
849 if (isHadron(parent)&&!isBeam(parent)) {
850 if (!hadron) hadron = parent; // assumes linear hadron parentage
851 return true;
852 }
853 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
854 }
855 return fromHad;
856 }
857
860
861 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
862 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
863 decltype(thePart->end_vertex()) pVert(nullptr);
864 if (EndVert != nullptr) {
865 do {
866 bool samePart = false;
867 pVert = nullptr;
868 auto outgoing = EndVert->particles_out();
869 auto incoming = EndVert->particles_in();
870 for (const auto& itrDaug: outgoing) {
871 if (!itrDaug) continue;
872 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
873 // brem on generator level for tau
874 (outgoing.size() == 1 && incoming.size() == 1 &&
876 itrDaug->pdg_id() == thePart->pdg_id()) {
877 samePart = true;
878 pVert = itrDaug->end_vertex();
879 }
880 }
881 if (samePart) EndVert = pVert;
882 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
883 }
884 return EndVert;
885 }
886
888
889 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
890 if (!theVert) return {};
891 decltype(theVert->particles_out()) finalStatePart;
892 auto outgoing = theVert->particles_out();
893 for (const auto& thePart: outgoing) {
894 if (!thePart) continue;
895 finalStatePart.push_back(thePart);
896 if (isStable(thePart)) continue;
897 V pVert = findSimulatedEndVertex(thePart);
898 if (pVert == theVert) break; // to prevent Sherpa loop
899 if (pVert != nullptr) {
900 auto vecPart = findFinalStateParticles<V>(pVert);
901 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
902 }
903 }
904 return finalStatePart;
905 }
906#if !defined(XAOD_ANALYSIS)
907#include "AtlasHepMC/GenEvent.h"
908#ifdef HEPMC3
909inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
910inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
911#else
912inline void GeVToMeV(HepMC::GenEvent* evt) {
913 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
914 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
915 (*p)->momentum().py() * 1000,
916 (*p)->momentum().pz() * 1000,
917 (*p)->momentum().e() * 1000);
918 (*p)->set_momentum(fv);
919 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
920 }
921}
922inline void MeVToGeV(HepMC::GenEvent* evt) {
923 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
924 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
925 (*p)->momentum().py() / 1000,
926 (*p)->momentum().pz() / 1000,
927 (*p)->momentum().e() / 1000);
928 (*p)->set_momentum(fv);
929 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
930 }
931}
932#endif
933#endif
934}
935#endif

◆ suggest_barcode() [1/2]

template<class T>
bool MC::HepMC::suggest_barcode ( T & p,
int i )

Definition at line 670 of file HepMCHelpers.h.

691{
692inline
693auto particles_in (const HepMC::GenVertex* p) {
694 return std::ranges::subrange (p->particles_in_const_begin(),
695 p->particles_in_const_end());
696}
697}
698#endif
699
700namespace MC
701{
702 template <class VTX>
703 auto particles_in (const VTX* p) { return p->particles_in(); }
704 template <class VTX>
705 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
706
707 namespace Pythia8
708 {
710 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
711
712 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
713
714 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
715 }
716
717#include "AtlasPID.h"
718
720 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
721
723 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
724
726 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
727
729 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
730
732 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
733
735 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
736
738 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
739
741 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
742
744 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
745
747 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
748
752 template <class T> inline bool isStableOrSimDecayed(const T& p) {
753 const auto vertex = p->end_vertex();
754 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
755 }
756
758 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
759
761 template <class T> inline bool isSpecialNonInteracting(const T& p) {
762 const int apid = std::abs(p->pdg_id());
763 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
764 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
765 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
766 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
767 return false;
768 }
769
771
772 template <class T> T findMother(T thePart) {
773 auto partOriVert = thePart->production_vertex();
774 if (!partOriVert) return nullptr;
775
776 long partPDG = thePart->pdg_id();
777 long MotherPDG(0);
778
779 auto MothOriVert = partOriVert;
780 MothOriVert = nullptr;
781 T theMoth(nullptr);
782
783 size_t itr = 0;
784 do {
785 if (itr != 0) partOriVert = MothOriVert;
786 for ( const auto& p : particles_in(partOriVert) ) {
787 theMoth = p;
788 if (!theMoth) continue;
789 MotherPDG = theMoth->pdg_id();
790 MothOriVert = theMoth->production_vertex();
791 if (MotherPDG == partPDG) break;
792 }
793 itr++;
794 if (itr > 100) {
795 break;
796 }
797 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
798 MothOriVert != partOriVert);
799 return theMoth;
800 }
801
803
804 template <class C, class T> T findMatching(C TruthContainer, T p) {
805 T ptrPart = nullptr;
806 if (!p) return ptrPart;
807 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
808 for (T truthParticle : *TruthContainer) {
809 if (HepMC::is_sim_descendant(p,truthParticle)) {
810 ptrPart = truthParticle;
811 break;
812 }
813 }
814 }
815 else {
816 for (T truthParticle : TruthContainer) {
817 if (HepMC::is_sim_descendant(p,truthParticle)) {
818 ptrPart = truthParticle;
819 break;
820 }
821 }
822 }
823 return ptrPart;
824 }
826
827 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
828 auto prodVtx = thePart->production_vertex();
829 if (!prodVtx) return;
830 for (const auto& theMother: prodVtx->particles_in()) {
831 if (!theMother) continue;
832 allancestors.insert(theMother);
833 findParticleAncestors(theMother, allancestors);
834 }
835 }
836
838
839 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
840 auto endVtx = thePart->end_vertex();
841 if (!endVtx) return;
842 for (const auto& theDaughter: endVtx->particles_out()) {
843 if (!theDaughter) continue;
844 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
845 allstabledescendants.insert(theDaughter);
846 }
847 findParticleStableDescendants(theDaughter, allstabledescendants);
848 }
849 }
850
854
855 template <class T> bool isHardScatteringVertex(T pVert) {
856 if (pVert == nullptr) return false;
857 T pV = pVert;
858 int numOfPartIn(0);
859 int pdg(0);
860
861 do {
862 pVert = pV;
863 auto incoming = pVert->particles_in();
864 numOfPartIn = incoming.size();
865 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
866 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
867
868 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
869
870 if (numOfPartIn == 2) {
871 auto incoming = pVert->particles_in();
872 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
873 }
874 return false;
875}
876
880
881 template <class T, class U>
882 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
883 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
884 auto vtx = p->production_vertex();
885 if (!vtx) return false;
886 bool fromHad = false;
887 for ( const auto& parent : particles_in(vtx) ) {
888 if (!parent) continue;
889 // should this really go into parton-level territory?
890 // probably depends where BSM particles are being decayed
891 fromBSM |= isBSM(parent);
892 if (!isPhysical(parent)) return false;
893 fromTau |= isTau(parent);
894 if (isHadron(parent)&&!isBeam(parent)) {
895 if (!hadron) hadron = parent; // assumes linear hadron parentage
896 return true;
897 }
898 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
899 }
900 return fromHad;
901 }
902
905
906 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
907 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
908 decltype(thePart->end_vertex()) pVert(nullptr);
909 if (EndVert != nullptr) {
910 do {
911 bool samePart = false;
912 pVert = nullptr;
913 auto outgoing = EndVert->particles_out();
914 auto incoming = EndVert->particles_in();
915 for (const auto& itrDaug: outgoing) {
916 if (!itrDaug) continue;
917 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
918 // brem on generator level for tau
919 (outgoing.size() == 1 && incoming.size() == 1 &&
921 itrDaug->pdg_id() == thePart->pdg_id()) {
922 samePart = true;
923 pVert = itrDaug->end_vertex();
924 }
925 }
926 if (samePart) EndVert = pVert;
927 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
928 }
929 return EndVert;
930 }
931
933
934 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
935 if (!theVert) return {};
936 decltype(theVert->particles_out()) finalStatePart;
937 auto outgoing = theVert->particles_out();
938 for (const auto& thePart: outgoing) {
939 if (!thePart) continue;
940 finalStatePart.push_back(thePart);
941 if (isStable(thePart)) continue;
942 V pVert = findSimulatedEndVertex(thePart);
943 if (pVert == theVert) break; // to prevent Sherpa loop
944 if (pVert != nullptr) {
945 auto vecPart = findFinalStateParticles<V>(pVert);
946 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
947 }
948 }
949 return finalStatePart;
950 }
951#if !defined(XAOD_ANALYSIS)
952#include "AtlasHepMC/GenEvent.h"
953#ifdef HEPMC3
954inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
955inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
956#else
957inline void GeVToMeV(HepMC::GenEvent* evt) {
958 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
959 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
960 (*p)->momentum().py() * 1000,
961 (*p)->momentum().pz() * 1000,
962 (*p)->momentum().e() * 1000);
963 (*p)->set_momentum(fv);
964 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
965 }
966}
967inline void MeVToGeV(HepMC::GenEvent* evt) {
968 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
969 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
970 (*p)->momentum().py() / 1000,
971 (*p)->momentum().pz() / 1000,
972 (*p)->momentum().e() / 1000);
973 (*p)->set_momentum(fv);
974 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
975 }
976}
977#endif
978#endif
979}
980#endif

◆ suggest_barcode() [2/2]

template<class T>
bool MC::HepMC::suggest_barcode ( T * p,
int i )

Definition at line 671 of file HepMCHelpers.h.

692{
693inline
694auto particles_in (const HepMC::GenVertex* p) {
695 return std::ranges::subrange (p->particles_in_const_begin(),
696 p->particles_in_const_end());
697}
698}
699#endif
700
701namespace MC
702{
703 template <class VTX>
704 auto particles_in (const VTX* p) { return p->particles_in(); }
705 template <class VTX>
706 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
707
708 namespace Pythia8
709 {
711 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
712
713 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
714
715 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
716 }
717
718#include "AtlasPID.h"
719
721 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
722
724 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
725
727 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
728
730 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
731
733 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
734
736 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
737
739 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
740
742 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
743
745 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
746
748 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
749
753 template <class T> inline bool isStableOrSimDecayed(const T& p) {
754 const auto vertex = p->end_vertex();
755 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
756 }
757
759 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
760
762 template <class T> inline bool isSpecialNonInteracting(const T& p) {
763 const int apid = std::abs(p->pdg_id());
764 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
765 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
766 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
767 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
768 return false;
769 }
770
772
773 template <class T> T findMother(T thePart) {
774 auto partOriVert = thePart->production_vertex();
775 if (!partOriVert) return nullptr;
776
777 long partPDG = thePart->pdg_id();
778 long MotherPDG(0);
779
780 auto MothOriVert = partOriVert;
781 MothOriVert = nullptr;
782 T theMoth(nullptr);
783
784 size_t itr = 0;
785 do {
786 if (itr != 0) partOriVert = MothOriVert;
787 for ( const auto& p : particles_in(partOriVert) ) {
788 theMoth = p;
789 if (!theMoth) continue;
790 MotherPDG = theMoth->pdg_id();
791 MothOriVert = theMoth->production_vertex();
792 if (MotherPDG == partPDG) break;
793 }
794 itr++;
795 if (itr > 100) {
796 break;
797 }
798 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
799 MothOriVert != partOriVert);
800 return theMoth;
801 }
802
804
805 template <class C, class T> T findMatching(C TruthContainer, T p) {
806 T ptrPart = nullptr;
807 if (!p) return ptrPart;
808 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
809 for (T truthParticle : *TruthContainer) {
810 if (HepMC::is_sim_descendant(p,truthParticle)) {
811 ptrPart = truthParticle;
812 break;
813 }
814 }
815 }
816 else {
817 for (T truthParticle : TruthContainer) {
818 if (HepMC::is_sim_descendant(p,truthParticle)) {
819 ptrPart = truthParticle;
820 break;
821 }
822 }
823 }
824 return ptrPart;
825 }
827
828 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
829 auto prodVtx = thePart->production_vertex();
830 if (!prodVtx) return;
831 for (const auto& theMother: prodVtx->particles_in()) {
832 if (!theMother) continue;
833 allancestors.insert(theMother);
834 findParticleAncestors(theMother, allancestors);
835 }
836 }
837
839
840 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
841 auto endVtx = thePart->end_vertex();
842 if (!endVtx) return;
843 for (const auto& theDaughter: endVtx->particles_out()) {
844 if (!theDaughter) continue;
845 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
846 allstabledescendants.insert(theDaughter);
847 }
848 findParticleStableDescendants(theDaughter, allstabledescendants);
849 }
850 }
851
855
856 template <class T> bool isHardScatteringVertex(T pVert) {
857 if (pVert == nullptr) return false;
858 T pV = pVert;
859 int numOfPartIn(0);
860 int pdg(0);
861
862 do {
863 pVert = pV;
864 auto incoming = pVert->particles_in();
865 numOfPartIn = incoming.size();
866 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
867 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
868
869 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
870
871 if (numOfPartIn == 2) {
872 auto incoming = pVert->particles_in();
873 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
874 }
875 return false;
876}
877
881
882 template <class T, class U>
883 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
884 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
885 auto vtx = p->production_vertex();
886 if (!vtx) return false;
887 bool fromHad = false;
888 for ( const auto& parent : particles_in(vtx) ) {
889 if (!parent) continue;
890 // should this really go into parton-level territory?
891 // probably depends where BSM particles are being decayed
892 fromBSM |= isBSM(parent);
893 if (!isPhysical(parent)) return false;
894 fromTau |= isTau(parent);
895 if (isHadron(parent)&&!isBeam(parent)) {
896 if (!hadron) hadron = parent; // assumes linear hadron parentage
897 return true;
898 }
899 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
900 }
901 return fromHad;
902 }
903
906
907 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
908 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
909 decltype(thePart->end_vertex()) pVert(nullptr);
910 if (EndVert != nullptr) {
911 do {
912 bool samePart = false;
913 pVert = nullptr;
914 auto outgoing = EndVert->particles_out();
915 auto incoming = EndVert->particles_in();
916 for (const auto& itrDaug: outgoing) {
917 if (!itrDaug) continue;
918 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
919 // brem on generator level for tau
920 (outgoing.size() == 1 && incoming.size() == 1 &&
922 itrDaug->pdg_id() == thePart->pdg_id()) {
923 samePart = true;
924 pVert = itrDaug->end_vertex();
925 }
926 }
927 if (samePart) EndVert = pVert;
928 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
929 }
930 return EndVert;
931 }
932
934
935 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
936 if (!theVert) return {};
937 decltype(theVert->particles_out()) finalStatePart;
938 auto outgoing = theVert->particles_out();
939 for (const auto& thePart: outgoing) {
940 if (!thePart) continue;
941 finalStatePart.push_back(thePart);
942 if (isStable(thePart)) continue;
943 V pVert = findSimulatedEndVertex(thePart);
944 if (pVert == theVert) break; // to prevent Sherpa loop
945 if (pVert != nullptr) {
946 auto vecPart = findFinalStateParticles<V>(pVert);
947 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
948 }
949 }
950 return finalStatePart;
951 }
952#if !defined(XAOD_ANALYSIS)
953#include "AtlasHepMC/GenEvent.h"
954#ifdef HEPMC3
955inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
956inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
957#else
958inline void GeVToMeV(HepMC::GenEvent* evt) {
959 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
960 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
961 (*p)->momentum().py() * 1000,
962 (*p)->momentum().pz() * 1000,
963 (*p)->momentum().e() * 1000);
964 (*p)->set_momentum(fv);
965 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
966 }
967}
968inline void MeVToGeV(HepMC::GenEvent* evt) {
969 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
970 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
971 (*p)->momentum().py() / 1000,
972 (*p)->momentum().pz() / 1000,
973 (*p)->momentum().e() / 1000);
974 (*p)->set_momentum(fv);
975 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
976 }
977}
978#endif
979#endif
980}
981#endif

◆ suggest_barcode< std::unique_ptr< HepMC::GenParticle > >()

template<>
bool MC::HepMC::suggest_barcode< std::unique_ptr< HepMC::GenParticle > > ( std::unique_ptr< HepMC::GenParticle > & p,
int i )
inline

Definition at line 673 of file HepMCHelpers.h.

694{
695inline
696auto particles_in (const HepMC::GenVertex* p) {
697 return std::ranges::subrange (p->particles_in_const_begin(),
698 p->particles_in_const_end());
699}
700}
701#endif
702
703namespace MC
704{
705 template <class VTX>
706 auto particles_in (const VTX* p) { return p->particles_in(); }
707 template <class VTX>
708 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
709
710 namespace Pythia8
711 {
713 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
714
715 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
716
717 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
718 }
719
720#include "AtlasPID.h"
721
723 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
724
726 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
727
729 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
730
732 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
733
735 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
736
738 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
739
741 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
742
744 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
745
747 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
748
750 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
751
755 template <class T> inline bool isStableOrSimDecayed(const T& p) {
756 const auto vertex = p->end_vertex();
757 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
758 }
759
761 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
762
764 template <class T> inline bool isSpecialNonInteracting(const T& p) {
765 const int apid = std::abs(p->pdg_id());
766 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
767 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
768 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
769 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
770 return false;
771 }
772
774
775 template <class T> T findMother(T thePart) {
776 auto partOriVert = thePart->production_vertex();
777 if (!partOriVert) return nullptr;
778
779 long partPDG = thePart->pdg_id();
780 long MotherPDG(0);
781
782 auto MothOriVert = partOriVert;
783 MothOriVert = nullptr;
784 T theMoth(nullptr);
785
786 size_t itr = 0;
787 do {
788 if (itr != 0) partOriVert = MothOriVert;
789 for ( const auto& p : particles_in(partOriVert) ) {
790 theMoth = p;
791 if (!theMoth) continue;
792 MotherPDG = theMoth->pdg_id();
793 MothOriVert = theMoth->production_vertex();
794 if (MotherPDG == partPDG) break;
795 }
796 itr++;
797 if (itr > 100) {
798 break;
799 }
800 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
801 MothOriVert != partOriVert);
802 return theMoth;
803 }
804
806
807 template <class C, class T> T findMatching(C TruthContainer, T p) {
808 T ptrPart = nullptr;
809 if (!p) return ptrPart;
810 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
811 for (T truthParticle : *TruthContainer) {
812 if (HepMC::is_sim_descendant(p,truthParticle)) {
813 ptrPart = truthParticle;
814 break;
815 }
816 }
817 }
818 else {
819 for (T truthParticle : TruthContainer) {
820 if (HepMC::is_sim_descendant(p,truthParticle)) {
821 ptrPart = truthParticle;
822 break;
823 }
824 }
825 }
826 return ptrPart;
827 }
829
830 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
831 auto prodVtx = thePart->production_vertex();
832 if (!prodVtx) return;
833 for (const auto& theMother: prodVtx->particles_in()) {
834 if (!theMother) continue;
835 allancestors.insert(theMother);
836 findParticleAncestors(theMother, allancestors);
837 }
838 }
839
841
842 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
843 auto endVtx = thePart->end_vertex();
844 if (!endVtx) return;
845 for (const auto& theDaughter: endVtx->particles_out()) {
846 if (!theDaughter) continue;
847 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
848 allstabledescendants.insert(theDaughter);
849 }
850 findParticleStableDescendants(theDaughter, allstabledescendants);
851 }
852 }
853
857
858 template <class T> bool isHardScatteringVertex(T pVert) {
859 if (pVert == nullptr) return false;
860 T pV = pVert;
861 int numOfPartIn(0);
862 int pdg(0);
863
864 do {
865 pVert = pV;
866 auto incoming = pVert->particles_in();
867 numOfPartIn = incoming.size();
868 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
869 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
870
871 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
872
873 if (numOfPartIn == 2) {
874 auto incoming = pVert->particles_in();
875 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
876 }
877 return false;
878}
879
883
884 template <class T, class U>
885 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
886 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
887 auto vtx = p->production_vertex();
888 if (!vtx) return false;
889 bool fromHad = false;
890 for ( const auto& parent : particles_in(vtx) ) {
891 if (!parent) continue;
892 // should this really go into parton-level territory?
893 // probably depends where BSM particles are being decayed
894 fromBSM |= isBSM(parent);
895 if (!isPhysical(parent)) return false;
896 fromTau |= isTau(parent);
897 if (isHadron(parent)&&!isBeam(parent)) {
898 if (!hadron) hadron = parent; // assumes linear hadron parentage
899 return true;
900 }
901 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
902 }
903 return fromHad;
904 }
905
908
909 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
910 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
911 decltype(thePart->end_vertex()) pVert(nullptr);
912 if (EndVert != nullptr) {
913 do {
914 bool samePart = false;
915 pVert = nullptr;
916 auto outgoing = EndVert->particles_out();
917 auto incoming = EndVert->particles_in();
918 for (const auto& itrDaug: outgoing) {
919 if (!itrDaug) continue;
920 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
921 // brem on generator level for tau
922 (outgoing.size() == 1 && incoming.size() == 1 &&
924 itrDaug->pdg_id() == thePart->pdg_id()) {
925 samePart = true;
926 pVert = itrDaug->end_vertex();
927 }
928 }
929 if (samePart) EndVert = pVert;
930 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
931 }
932 return EndVert;
933 }
934
936
937 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
938 if (!theVert) return {};
939 decltype(theVert->particles_out()) finalStatePart;
940 auto outgoing = theVert->particles_out();
941 for (const auto& thePart: outgoing) {
942 if (!thePart) continue;
943 finalStatePart.push_back(thePart);
944 if (isStable(thePart)) continue;
945 V pVert = findSimulatedEndVertex(thePart);
946 if (pVert == theVert) break; // to prevent Sherpa loop
947 if (pVert != nullptr) {
948 auto vecPart = findFinalStateParticles<V>(pVert);
949 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
950 }
951 }
952 return finalStatePart;
953 }
954#if !defined(XAOD_ANALYSIS)
955#include "AtlasHepMC/GenEvent.h"
956#ifdef HEPMC3
957inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
958inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
959#else
960inline void GeVToMeV(HepMC::GenEvent* evt) {
961 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
962 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
963 (*p)->momentum().py() * 1000,
964 (*p)->momentum().pz() * 1000,
965 (*p)->momentum().e() * 1000);
966 (*p)->set_momentum(fv);
967 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
968 }
969}
970inline void MeVToGeV(HepMC::GenEvent* evt) {
971 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
972 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
973 (*p)->momentum().py() / 1000,
974 (*p)->momentum().pz() / 1000,
975 (*p)->momentum().e() / 1000);
976 (*p)->set_momentum(fv);
977 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
978 }
979}
980#endif
981#endif
982}
983#endif

◆ valid_beam_particles()

bool MC::HepMC::valid_beam_particles ( const GenEvent * e)
inline

Definition at line 681 of file HepMCHelpers.h.

702{
703inline
704auto particles_in (const HepMC::GenVertex* p) {
705 return std::ranges::subrange (p->particles_in_const_begin(),
706 p->particles_in_const_end());
707}
708}
709#endif
710
711namespace MC
712{
713 template <class VTX>
714 auto particles_in (const VTX* p) { return p->particles_in(); }
715 template <class VTX>
716 auto particles_in (const std::shared_ptr<VTX>& p) { return p->particles_in(); }
717
718 namespace Pythia8
719 {
721 template <class T> inline bool isConditionA(const T& p) { return p->status() == 62 || p->status() == 52 || p->status() == 21 || p->status() == 22;}
722
723 template <class T> inline bool isConditionB(const T& p) { return p->status() == 23;}
724
725 template <class T> inline bool isConditionC(const T& p) { return p->status() > 30 && p->status() < 40;}
726 }
727
728#include "AtlasPID.h"
729
731 template <class T> inline bool isInteracting(const T& p) { return isStrongInteracting<T>(p) || isEMInteracting<T>(p) || isGeantino<T>(p); }
732
734 template <class T> inline bool isChargedNonShowering(const T& p) { return (isMuon<T>(p) || isSUSY<T>(p)); }
735
737 template <class T> inline bool isBeam(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 4;}
738
740 template <class T> inline bool isDecayed(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 2;}
741
743 template <class T> inline bool isStable(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1;}
744
746 template <class T> inline bool isFinalState(const T& p) { return HepMC::status(p)%HepMC::SIM_STATUS_THRESHOLD == 1 && !p->end_vertex();}
747
749 template <class T> inline bool isPhysical(const T& p) { return isStable<T>(p) || isDecayed<T>(p); }
750
752 template <class T> inline bool isGenStable(const T& p) { return isStable<T>(p) && !HepMC::is_simulation_particle<T>(p);}
753
755 template <class T> inline bool isSimStable(const T& p) { return isStable<T>(p) && !p->end_vertex() && HepMC::is_simulation_particle<T>(p);}
756
758 template <class T> inline bool isSimInteracting(const T& p) { return isGenStable<T>(p) && isInteracting<T>(p);}
759
763 template <class T> inline bool isStableOrSimDecayed(const T& p) {
764 const auto vertex = p->end_vertex();
765 return ( isStable<T>(p) || (isDecayed<T>(p) && (!vertex || HepMC::is_simulation_vertex(vertex))));
766 }
767
769 template <class T> inline bool isZeroEnergyPhoton(const T& p) { return isPhoton<T>(p) && p->e() == 0;}
770
772 template <class T> inline bool isSpecialNonInteracting(const T& p) {
773 const int apid = std::abs(p->pdg_id());
774 if (apid == NU_E || apid == NU_MU || apid == NU_TAU) return true; //< neutrinos
775 if (apid == 1000022 || apid == 1000024 || apid == 5100022) return true; // SUSY & KK photon and Z partners
776 if (apid == GRAVITON || apid == 1000039 || apid == 5000039) return true; //< gravitons: standard, SUSY and KK
777 if (apid == 9000001 || apid == 9000002 || apid == 9000003 || apid == 9000004 || apid == 9000005 || apid == 9000006) return true; //< exotic particles from monotop model
778 return false;
779 }
780
782
783 template <class T> T findMother(T thePart) {
784 auto partOriVert = thePart->production_vertex();
785 if (!partOriVert) return nullptr;
786
787 long partPDG = thePart->pdg_id();
788 long MotherPDG(0);
789
790 auto MothOriVert = partOriVert;
791 MothOriVert = nullptr;
792 T theMoth(nullptr);
793
794 size_t itr = 0;
795 do {
796 if (itr != 0) partOriVert = MothOriVert;
797 for ( const auto& p : particles_in(partOriVert) ) {
798 theMoth = p;
799 if (!theMoth) continue;
800 MotherPDG = theMoth->pdg_id();
801 MothOriVert = theMoth->production_vertex();
802 if (MotherPDG == partPDG) break;
803 }
804 itr++;
805 if (itr > 100) {
806 break;
807 }
808 } while (MothOriVert != nullptr && MotherPDG == partPDG && !HepMC::is_simulation_particle(thePart) &&
809 MothOriVert != partOriVert);
810 return theMoth;
811 }
812
814
815 template <class C, class T> T findMatching(C TruthContainer, T p) {
816 T ptrPart = nullptr;
817 if (!p) return ptrPart;
818 if constexpr (std::is_pointer_v<C> || HepMC::is_smart_ptr_v<C>){ //C is ptr
819 for (T truthParticle : *TruthContainer) {
820 if (HepMC::is_sim_descendant(p,truthParticle)) {
821 ptrPart = truthParticle;
822 break;
823 }
824 }
825 }
826 else {
827 for (T truthParticle : TruthContainer) {
828 if (HepMC::is_sim_descendant(p,truthParticle)) {
829 ptrPart = truthParticle;
830 break;
831 }
832 }
833 }
834 return ptrPart;
835 }
837
838 template <class T> void findParticleAncestors(T thePart, std::set<T>& allancestors) {
839 auto prodVtx = thePart->production_vertex();
840 if (!prodVtx) return;
841 for (const auto& theMother: prodVtx->particles_in()) {
842 if (!theMother) continue;
843 allancestors.insert(theMother);
844 findParticleAncestors(theMother, allancestors);
845 }
846 }
847
849
850 template <class T> void findParticleStableDescendants(T thePart, std::set<T>& allstabledescendants) {
851 auto endVtx = thePart->end_vertex();
852 if (!endVtx) return;
853 for (const auto& theDaughter: endVtx->particles_out()) {
854 if (!theDaughter) continue;
855 if (isStable(theDaughter) && !HepMC::is_simulation_particle(theDaughter)) {
856 allstabledescendants.insert(theDaughter);
857 }
858 findParticleStableDescendants(theDaughter, allstabledescendants);
859 }
860 }
861
865
866 template <class T> bool isHardScatteringVertex(T pVert) {
867 if (pVert == nullptr) return false;
868 T pV = pVert;
869 int numOfPartIn(0);
870 int pdg(0);
871
872 do {
873 pVert = pV;
874 auto incoming = pVert->particles_in();
875 numOfPartIn = incoming.size();
876 pdg = numOfPartIn && incoming.front() != nullptr ? incoming.front()->pdg_id() : 0;
877 pV = numOfPartIn && incoming.front() != nullptr ? incoming.front()->production_vertex() : nullptr;
878
879 } while (numOfPartIn == 1 && (std::abs(pdg) < 81 || std::abs(pdg) > 100) && pV != nullptr);
880
881 if (numOfPartIn == 2) {
882 auto incoming = pVert->particles_in();
883 if (incoming.at(0) && incoming.at(1) && (std::abs(incoming.at(0)->pdg_id()) < 7 || incoming.at(0)->pdg_id() == 21) && (std::abs(incoming.at(1)->pdg_id()) < 7 || incoming.at(1)->pdg_id() == 21)) return true;
884 }
885 return false;
886}
887
891
892 template <class T, class U>
893 bool isFromHadron(T p, U hadron, bool &fromTau, bool &fromBSM) {
894 if (isHadron(p)&&!isBeam(p)) return true; // trivial case
895 auto vtx = p->production_vertex();
896 if (!vtx) return false;
897 bool fromHad = false;
898 for ( const auto& parent : particles_in(vtx) ) {
899 if (!parent) continue;
900 // should this really go into parton-level territory?
901 // probably depends where BSM particles are being decayed
902 fromBSM |= isBSM(parent);
903 if (!isPhysical(parent)) return false;
904 fromTau |= isTau(parent);
905 if (isHadron(parent)&&!isBeam(parent)) {
906 if (!hadron) hadron = parent; // assumes linear hadron parentage
907 return true;
908 }
909 fromHad |= isFromHadron(parent, hadron, fromTau, fromBSM);
910 }
911 return fromHad;
912 }
913
916
917 template <class T> auto findSimulatedEndVertex(T thePart) -> decltype(thePart->end_vertex()) {
918 decltype(thePart->end_vertex()) EndVert = thePart->end_vertex();
919 decltype(thePart->end_vertex()) pVert(nullptr);
920 if (EndVert != nullptr) {
921 do {
922 bool samePart = false;
923 pVert = nullptr;
924 auto outgoing = EndVert->particles_out();
925 auto incoming = EndVert->particles_in();
926 for (const auto& itrDaug: outgoing) {
927 if (!itrDaug) continue;
928 if ((( HepMC::is_same_generator_particle(itrDaug,thePart)) ||
929 // brem on generator level for tau
930 (outgoing.size() == 1 && incoming.size() == 1 &&
932 itrDaug->pdg_id() == thePart->pdg_id()) {
933 samePart = true;
934 pVert = itrDaug->end_vertex();
935 }
936 }
937 if (samePart) EndVert = pVert;
938 } while (pVert != nullptr && pVert != EndVert); // pVert!=EndVert to prevent Sherpa loop
939 }
940 return EndVert;
941 }
942
944
945 template <class V> auto findFinalStateParticles(V theVert) -> decltype(theVert->particles_out()) {
946 if (!theVert) return {};
947 decltype(theVert->particles_out()) finalStatePart;
948 auto outgoing = theVert->particles_out();
949 for (const auto& thePart: outgoing) {
950 if (!thePart) continue;
951 finalStatePart.push_back(thePart);
952 if (isStable(thePart)) continue;
953 V pVert = findSimulatedEndVertex(thePart);
954 if (pVert == theVert) break; // to prevent Sherpa loop
955 if (pVert != nullptr) {
956 auto vecPart = findFinalStateParticles<V>(pVert);
957 finalStatePart.insert(finalStatePart.end(),vecPart.begin(),vecPart.end());
958 }
959 }
960 return finalStatePart;
961 }
962#if !defined(XAOD_ANALYSIS)
963#include "AtlasHepMC/GenEvent.h"
964#ifdef HEPMC3
965inline void GeVToMeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1000); p->set_generated_mass(1000* p->generated_mass());}}
966inline void MeVToGeV(HepMC::GenEvent* evt) { for (auto& p: evt->particles()) { p->set_momentum(p->momentum()*1.0/1000); p->set_generated_mass(1.0/1000* p->generated_mass());} }
967#else
968inline void GeVToMeV(HepMC::GenEvent* evt) {
969 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
970 const HepMC::FourVector fv((*p)->momentum().px() * 1000,
971 (*p)->momentum().py() * 1000,
972 (*p)->momentum().pz() * 1000,
973 (*p)->momentum().e() * 1000);
974 (*p)->set_momentum(fv);
975 (*p)->set_generated_mass(1000 * (*p)->generated_mass());
976 }
977}
978inline void MeVToGeV(HepMC::GenEvent* evt) {
979 for (HepMC::GenEvent::particle_iterator p = evt->particles_begin(); p != evt->particles_end(); ++p) {
980 const HepMC::FourVector fv((*p)->momentum().px() / 1000,
981 (*p)->momentum().py() / 1000,
982 (*p)->momentum().pz() / 1000,
983 (*p)->momentum().e() / 1000);
984 (*p)->set_momentum(fv);
985 (*p)->set_generated_mass((*p)->generated_mass() / 1000);
986 }
987}
988#endif
989#endif
990}
991#endif