2 #ifndef RIVET_Particle_HH 3 #define RIVET_Particle_HH 5 #include "Rivet/Particle.fhh" 6 #include "Rivet/ParticleBase.hh" 7 #include "Rivet/Config/RivetCommon.hh" 8 #include "Rivet/Tools/Cuts.hh" 9 #include "Rivet/Tools/Utils.hh" 10 #include "Rivet/Math/LorentzTrans.hh" 12 #include "fastjet/PseudoJet.hh" 28 _original(nullptr), _id(PID::ANY), _isDirect(4,
std::make_pair(false,false))
34 _original(gp), _id(pid),
35 _momentum(mom), _origin(pos),
36 _isDirect(4,
std::make_pair(false,false))
47 _original(gp), _id(gp->pdg_id()),
49 _isDirect(4,
std::make_pair(false,false))
51 ConstGenVertexPtr vprod = gp->production_vertex();
52 if (vprod !=
nullptr) {
53 setOrigin(vprod->position().t(), vprod->position().x(), vprod->position().y(), vprod->position().z());
59 :
Particle(HepMCUtils::getParticlePtr(gp))
118 const double phi0 = M_PI/2 - this->
phi();
119 const double x0 = rho0 * cos(phi0);
120 const double y0 = rho0 * sin(phi0);
153 explicit operator ConstGenParticlePtr ()
const {
return genParticle(); }
162 PdgId
pid()
const {
return _id; }
164 PdgId
abspid()
const {
return std::abs(_id); }
266 const Particles
constituents(
const Cut& c,
const ParticleSorter& sorter)
const {
278 const Particles
constituents(
const ParticleSelector& selector,
const ParticleSorter& sorter)
const {
313 const Particles
rawConstituents(
const ParticleSelector& selector,
const ParticleSorter& sorter)
const {
328 Particles
parents(
const Cut& c=Cuts::OPEN)
const;
335 Particles
parents(
const ParticleSelector& f)
const {
386 Particles
ancestors(
const Cut& c=Cuts::OPEN,
bool only_physical=
true)
const;
394 Particles
ancestors(
const ParticleSelector& f,
bool only_physical=
true)
const {
404 return !
ancestors(f, only_physical).empty();
435 bool hasAncestor(PdgId pid,
bool only_physical=
true)
const;
477 bool fromTau(
bool prompt_taus_only=
false)
const;
502 DEPRECATED(
"Too vague: use fromHadron() || fromPromptTau(), or isDirect()")
515 bool isDirect(
bool allow_from_direct_tau=
false,
bool allow_from_direct_mu=
false)
const;
518 bool isPrompt(
bool allow_from_prompt_tau=
false,
bool allow_from_prompt_mu=
false)
const {
519 return isDirect(allow_from_prompt_tau, allow_from_prompt_mu);
535 Particles
children(
const Cut& c=Cuts::OPEN)
const;
538 Particles
children(
const ParticleSelector& f)
const {
574 Particles
allDescendants(
const Cut& c=Cuts::OPEN,
bool remove_duplicates=
true)
const;
577 Particles
allDescendants(
const ParticleSelector& f,
bool remove_duplicates=
true)
const {
665 if (!f(*
this))
return false;
677 if (!f(*
this))
return false;
697 if (
pid() != other.
pid())
return false;
709 ConstGenParticlePtr _original;
712 Particles _constituents;
725 mutable std::vector<std::pair<bool,bool> > _isDirect;
737 std::ostream&
operator << (std::ostream& os,
const ParticlePair& pp);
744 #include "Rivet/Tools/ParticleUtils.hh" virtual void setConstituents(const Particles &cs, bool setmom=false)
Set direct constituents of this particle.
Definition: MC_Cent_pPb.hh:10
int abscharge3(int pid)
Return the absolute value of 3 times the EM charge.
Definition: ParticleIdUtils.hh:865
bool any(const CONTAINER &c)
Return true if x is true for any x in container c, otherwise false.
Definition: Utils.hh:323
const Particles constituents(const ParticleSelector &selector) const
Direct constituents of this particle, filtered by a selection functor.
Definition: Particle.hh:272
double perp() const
Synonym for polarRadius.
Definition: Vector4.hh:108
double flightLength() const
Flight length (divide by mm or cm to get the appropriate units)
double pz() const
z component of momentum.
Definition: ParticleBase.hh:124
Particles children(const ParticleSelector &f) const
Get a list of the direct descendants from the current particle (with selector function) ...
Definition: Particle.hh:538
bool fromHadron() const
Determine whether the particle is from a hadron decay.
Particles rawConstituents() const
Fundamental constituents of this particle.
bool hasParent(PdgId pid) const
bool isParton(int pid)
Determine if the PID is that of a parton (quark or gluon)
Definition: ParticleIdUtils.hh:132
Particles ancestors(const Cut &c=Cuts::OPEN, bool only_physical=true) const
bool fromCharm() const
Determine whether the particle is from a c-hadron decay.
const FourMomentum & mom() const
Get equivalent single momentum four-vector (const) (alias).
Definition: ParticleBase.hh:39
const Particles rawConstituents(const ParticleSorter &sorter) const
Fundamental constituents of this particle, sorted by a functor.
Definition: Particle.hh:289
Vector3 closestApproach() const
Find the point of closest approach to the primary vertex.
Definition: Particle.hh:114
bool isFirstWithout(const ParticleSelector &f) const
Determine whether a particle is the first in a decay chain not to meet the function requirement...
Definition: Particle.hh:671
Jets filter_select(const Jets &jets, const Cut &c)
Filter a jet collection in-place to the subset that passes the supplied Cut.
Definition: JetUtils.hh:160
const Particles constituents(const Cut &c, const ParticleSorter &sorter) const
Direct constituents of this particle, sorted by a functor.
Definition: Particle.hh:266
bool isParton() const
Is this a parton? (Hopefully not very often... fiducial FTW)
Definition: Particle.hh:224
bool hasDescendantWith(const ParticleSelector &f, bool remove_duplicates=true) const
Definition: Particle.hh:586
Particle(const RivetHepMC::GenParticle &gp)
Constructor from a HepMC GenParticle reference.
Definition: Particle.hh:58
bool isStable() const
Whether this particle is stable according to the generator.
bool isBaryon(int pid)
Check to see if this is a valid baryon.
Definition: ParticleIdUtils.hh:260
double charge(int pid)
Return the EM charge (as floating point)
Definition: ParticleIdUtils.hh:868
PdgId pid() const
This Particle's PDG ID code.
Definition: Particle.hh:162
bool hasAncestorWith(const ParticleSelector &f, bool only_physical=true) const
Definition: Particle.hh:403
DEPRECATED("Too vague: use fromHadron() || fromPromptTau(), or isDirect()") bool fromDecay() const
Determine whether the particle is from a hadron or tau decay.
Definition: Particle.hh:502
Particles children(const Cut &c=Cuts::OPEN) const
Get a list of the direct descendants from the current particle (with optional selection Cut) ...
Particle representation, either from a HepMC::GenEvent or reconstructed.
Definition: Particle.hh:18
const Particles constituents(const Cut &c) const
Direct constituents of this particle, filtered by a Cut.
Definition: Particle.hh:260
virtual void addConstituent(const Particle &c, bool addmom=false)
Add a single direct constituent to this particle.
bool isHadron() const
Is this a hadron?
Definition: Particle.hh:194
Particle & setMomentum(double E, double px, double py, double pz)
Set the momentum via components.
Definition: Particle.hh:79
Particles parents(const ParticleSelector &f) const
Definition: Particle.hh:335
bool isCharged() const
Is this Particle charged?
Definition: Particle.hh:185
double theta() const
Synonym for polarAngle.
Definition: ParticleBase.hh:136
bool hasStableDescendantWith(const ParticleSelector &f) const
Definition: Particle.hh:628
Particle & setOrigin(double t, double x, double y, double z)
Set the origin position via components.
Definition: Particle.hh:103
bool hasChildWith(const ParticleSelector &f) const
Definition: Particle.hh:547
int abscharge3() const
Three times the absolute charge of this Particle (i.e. integer multiple of smallest quark charge)...
Definition: Particle.hh:182
bool isHadron(int pid)
Definition: ParticleIdUtils.hh:322
double phi(const PhiMapping mapping=ZERO_2PI) const
Get the directly.
Definition: ParticleBase.hh:105
double E() const
Get the energy directly.
Definition: ParticleBase.hh:51
const Particles & constituents() const
Direct constituents of this particle, returned by reference.
Definition: Particle.hh:250
Particles ancestors(const ParticleSelector &f, bool only_physical=true) const
Definition: Particle.hh:394
Particles stableDescendants(const ParticleSelector &f) const
Get a list of all the stable descendants from the current particle (with selector function) ...
Definition: Particle.hh:619
bool hasStableDescendantWithout(const ParticleSelector &f) const
Definition: Particle.hh:643
int charge3(int pid)
Three times the EM charge (as integer)
Definition: ParticleIdUtils.hh:800
bool isChargedLepton() const
Is this a charged lepton?
Definition: Particle.hh:206
virtual void addConstituents(const Particles &cs, bool addmom=false)
Add direct constituents to this particle.
double charge() const
The charge of this Particle.
Definition: Particle.hh:173
bool fromPromptTau() const
Determine whether the particle is from a prompt tau decay.
Definition: Particle.hh:484
Particle(PdgId pid, const FourMomentum &mom, const FourVector &pos=FourVector(), ConstGenParticlePtr gp=nullptr)
Constructor from PID and momentum.
Definition: Particle.hh:32
bool isLepton(int pid)
Definition: ParticleIdUtils.hh:340
MOMS sortBy(const MOMS &pbs, const CMP &cmp)
Sort a container of momenta by cmp and return by value for const inputs.
Definition: Vector4.hh:1437
double abscharge(int pid)
Return the EM charge (as floating point)
Definition: ParticleIdUtils.hh:871
const Particles constituents(const ParticleSorter &sorter) const
Direct constituents of this particle, sorted by a functor.
Definition: Particle.hh:254
const Particles rawConstituents(const Cut &c) const
Fundamental constituents of this particle, filtered by a Cut.
Definition: Particle.hh:295
const FourVector & origin() const
The origin position (and time).
Definition: Particle.hh:94
bool hasBottom(int pid)
Does this particle contain a bottom quark?
Definition: ParticleIdUtils.hh:581
Particle(ConstGenParticlePtr gp)
Constructor from a HepMC GenParticle pointer.
Definition: Particle.hh:45
bool isBaryon() const
Is this a baryon?
Definition: Particle.hh:200
bool isVisible() const
Is this particle potentially visible in a detector?
bool fromHadronicTau(bool prompt_taus_only=false) const
Determine whether the particle is from a tau which decayed hadronically.
bool isSame(const Particle &other) const
Definition: Particle.hh:696
bool hasCharm(int pid)
Does this particle contain a charm quark?
Definition: ParticleIdUtils.hh:579
Specialisation of VectorN to a general (non-momentum) Lorentz 4-vector.
Definition: Vector4.hh:22
bool isNeutrino() const
Is this a neutrino?
Definition: Particle.hh:209
const Particles rawConstituents(const ParticleSelector &selector, const ParticleSorter &sorter) const
Fundamental constituents of this particle, filtered and sorted by functors.
Definition: Particle.hh:313
const FourMomentum & momentum() const
The momentum.
Definition: Particle.hh:68
Particles allDescendants(const ParticleSelector &f, bool remove_duplicates=true) const
Get a list of all the descendants from the current particle (with selector function) ...
Definition: Particle.hh:577
Base class for particle-like things like Particle and Jet.
Definition: ParticleBase.hh:13
bool isLastWith(const ParticleSelector &f) const
Determine whether a particle is the last in a decay chain to meet the function requirement.
Definition: Particle.hh:676
bool hasDescendantWithout(const ParticleSelector &f, bool remove_duplicates=true) const
Definition: Particle.hh:601
bool isMeson() const
Is this a meson?
Definition: Particle.hh:197
bool isComposite() const
Determine if this Particle is a composite of other Rivet Particles.
Definition: Particle.hh:243
int charge3() const
Three times the charge of this Particle (i.e. integer multiple of smallest quark charge).
Definition: Particle.hh:179
Particle()
Definition: Particle.hh:26
Particle(PdgId pid, const FourMomentum &mom, ConstGenParticlePtr gp, const FourVector &pos=FourVector())
Constructor from PID, momentum, and a GenParticle for relational links.
Definition: Particle.hh:40
bool hasAncestorWithout(const ParticleSelector &f, bool only_physical=true) const
Definition: Particle.hh:418
bool isLastWithout(const ParticleSelector &f) const
Determine whether a particle is the last in a decay chain not to meet the function requirement...
Definition: Particle.hh:683
Particle & transformBy(const LorentzTransform <)
Apply an active Lorentz transform to this particle.
double px() const
x component of momentum.
Definition: ParticleBase.hh:120
bool isLepton() const
Is this a lepton?
Definition: Particle.hh:203
const Particles rawConstituents(const ParticleSelector &selector) const
Fundamental constituents of this particle, filtered by a selection functor.
Definition: Particle.hh:307
std::enable_if< std::is_floating_point< NUM >::value, bool >::type isZero(NUM val, double tolerance=1e-8)
Compare a number to zero.
Definition: MathUtils.hh:24
const Particles rawConstituents(const Cut &c, const ParticleSorter &sorter) const
Fundamental constituents of this particle, sorted by a functor.
Definition: Particle.hh:301
bool isFirstWith(const ParticleSelector &f) const
Determine whether a particle is the first in a decay chain to meet the function requirement.
Definition: Particle.hh:664
Particle & setMomentum(const FourMomentum &momentum)
Set the momentum.
Definition: Particle.hh:73
ConstGenParticlePtr genParticle() const
Get a const pointer to the original GenParticle.
Definition: Particle.hh:147
bool hasChildWithout(const ParticleSelector &f) const
Definition: Particle.hh:562
double abscharge() const
The absolute charge of this Particle.
Definition: Particle.hh:176
PdgId abspid() const
Absolute value of the PDG ID code.
Definition: Particle.hh:164
bool isMeson(int pid)
Check to see if this is a valid meson.
Definition: ParticleIdUtils.hh:237
bool fromTau(bool prompt_taus_only=false) const
Determine whether the particle is from a tau decay.
Three-dimensional specialisation of Vector.
Definition: Vector3.hh:26
bool hasBottom() const
Does this (hadron) contain a b quark?
Definition: Particle.hh:212
bool isDirect(bool allow_from_direct_tau=false, bool allow_from_direct_mu=false) const
Shorthand definition of 'promptness' based on set definition flags.
bool isNeutrino(int pid)
Determine if the PID is that of a neutrino.
Definition: ParticleIdUtils.hh:167
bool hasParentWithout(const ParticleSelector &f) const
Definition: Particle.hh:359
std::ostream & operator<<(std::ostream &os, const AnalysisInfo &ai)
Stream an AnalysisInfo as a text description.
Definition: AnalysisInfo.hh:365
Particles stableDescendants(const Cut &c=Cuts::OPEN) const
const Particles constituents(const ParticleSelector &selector, const ParticleSorter &sorter) const
Direct constituents of this particle, filtered and sorted by functors.
Definition: Particle.hh:278
Specialized version of the FourVector with momentum/energy functionality.
Definition: Vector4.hh:301
bool hasCharm() const
Does this (hadron) contain a c quark?
Definition: Particle.hh:215
double p() const
Get the 3-momentum magnitude directly.
Definition: ParticleBase.hh:110
bool isPrompt(bool allow_from_prompt_tau=false, bool allow_from_prompt_mu=false) const
Alias for isDirect.
Definition: Particle.hh:518
Particles parents(const Cut &c=Cuts::OPEN) const
Particle & setOrigin(const FourVector &position)
Set the origin position.
Definition: Particle.hh:98
double py() const
y component of momentum.
Definition: ParticleBase.hh:122
bool hasAncestor(PdgId pid, bool only_physical=true) const
virtual fastjet::PseudoJet pseudojet() const
Converter to FastJet3 PseudoJet.
Definition: Particle.hh:132
bool hasParentWith(const ParticleSelector &f) const
Definition: Particle.hh:344
Particle & setGenParticle(ConstGenParticlePtr gp)
Set a const pointer to the original GenParticle.
Definition: Particle.hh:141
bool isChargedLepton(int pid)
Determine if the PID is that of a charged lepton.
Definition: ParticleIdUtils.hh:158
Particles allDescendants(const Cut &c=Cuts::OPEN, bool remove_duplicates=true) const
Get a list of all the descendants from the current particle (with optional selection Cut) ...
bool fromBottom() const
Determine whether the particle is from a b-hadron decay.