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 648 of file HepMCHelpers.h.

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

◆ begin() [1/3]

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

Definition at line 642 of file HepMCHelpers.h.

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

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

◆ copyemptyGenEvent()

GenEvent * MC::HepMC::copyemptyGenEvent ( const GenEvent * inEvt)
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

◆ end() [1/3]

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

Definition at line 643 of file HepMCHelpers.h.

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

◆ fillBarcodesAttribute()

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

Definition at line 646 of file HepMCHelpers.h.

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

◆ get_ll_event_number()

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

Definition at line 637 of file HepMCHelpers.h.

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

◆ mpi() [1/2]

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

Definition at line 649 of file HepMCHelpers.h.

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

◆ mpi() [2/2]

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

Definition at line 652 of file HepMCHelpers.h.

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

◆ newGenEvent()

GenEvent * MC::HepMC::newGenEvent ( const int a,
const int b )
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

◆ 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 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

◆ set_mpi()

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

Definition at line 664 of file HepMCHelpers.h.

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

◆ set_random_states()

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

Definition at line 667 of file HepMCHelpers.h.

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

◆ set_signal_process_id()

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

Definition at line 661 of file HepMCHelpers.h.

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

◆ set_signal_process_vertex()

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

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

◆ signal_process_id() [1/2]

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

Definition at line 655 of file HepMCHelpers.h.

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

◆ signal_process_id() [2/2]

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

Definition at line 658 of file HepMCHelpers.h.

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

◆ signal_process_vertex()

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

Definition at line 645 of file HepMCHelpers.h.

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

◆ suggest_barcode() [1/2]

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

Definition at line 690 of file HepMCHelpers.h.

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

◆ suggest_barcode() [2/2]

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

Definition at line 691 of file HepMCHelpers.h.

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

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

◆ valid_beam_particles()

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

Definition at line 701 of file HepMCHelpers.h.

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