2 #ifndef RIVET_Particle_HH
3 #define RIVET_Particle_HH
5 #include "Rivet/Rivet.hh"
6 #include "Rivet/Particle.fhh"
7 #include "Rivet/ParticleBase.hh"
8 #include "Rivet/ParticleName.hh"
9 #include "Rivet/Math/Vectors.hh"
10 #include "Rivet/Tools/Logging.hh"
23 _original(0), _id(0), _momentum()
29 _original(0), _id(pid), _momentum(mom)
35 _original(&gp), _id(gp.pdg_id()),
51 return bool(_original);
118 const GenParticle* _original;
136 << pair.first.momentum().E()/GeV <<
" GeV, "
138 << pair.second.momentum().E()/GeV <<
" GeV]";
155 return a.momentum().pT() > b.momentum().pT();
219 inline double deltaR(
const Particle& p1,
const Particle& p2,
221 return deltaR(p1.momentum(), p2.momentum(), scheme);
224 inline double deltaR(
const Particle& p,
const FourMomentum& v,
226 return deltaR(p.momentum(), v, scheme);
229 inline double deltaR(
const Particle& p,
const FourVector& v,
231 return deltaR(p.momentum(), v, scheme);
234 inline double deltaR(
const Particle& p,
const Vector3& v) {
235 return deltaR(p.momentum(), v);
238 inline double deltaR(
const Particle& p,
double eta,
double phi) {
239 return deltaR(p.momentum(),
eta,
phi);
242 inline double deltaR(
const FourMomentum& v,
const Particle& p,
244 return deltaR(v, p.momentum(), scheme);
247 inline double deltaR(
const FourVector& v,
const Particle& p,
249 return deltaR(v, p.momentum(), scheme);
252 inline double deltaR(
const Vector3& v,
const Particle& p) {
253 return deltaR(v, p.momentum());
256 inline double deltaR(
double eta,
double phi,
const Particle& p) {
257 return deltaR(eta, phi, p.momentum());
261 inline double deltaPhi(
const Particle& p1,
const Particle& p2) {
262 return deltaPhi(p1.momentum(), p2.momentum());
265 inline double deltaPhi(
const Particle& p,
const FourMomentum& v) {
266 return deltaPhi(p.momentum(), v);
269 inline double deltaPhi(
const Particle& p,
const FourVector& v) {
270 return deltaPhi(p.momentum(), v);
273 inline double deltaPhi(
const Particle& p,
const Vector3& v) {
274 return deltaPhi(p.momentum(), v);
277 inline double deltaPhi(
const Particle& p,
double phi) {
278 return deltaPhi(p.momentum(),
phi);
281 inline double deltaPhi(
const FourMomentum& v,
const Particle& p) {
282 return deltaPhi(v, p.momentum());
285 inline double deltaPhi(
const FourVector& v,
const Particle& p) {
286 return deltaPhi(v, p.momentum());
289 inline double deltaPhi(
const Vector3& v,
const Particle& p) {
290 return deltaPhi(v, p.momentum());
293 inline double deltaPhi(
double phi,
const Particle& p) {
294 return deltaPhi(phi, p.momentum());
298 inline double deltaEta(
const Particle& p1,
const Particle& p2) {
299 return deltaEta(p1.momentum(), p2.momentum());
302 inline double deltaEta(
const Particle& p,
const FourMomentum& v) {
303 return deltaEta(p.momentum(), v);
306 inline double deltaEta(
const Particle& p,
const FourVector& v) {
307 return deltaEta(p.momentum(), v);
310 inline double deltaEta(
const Particle& p,
const Vector3& v) {
311 return deltaEta(p.momentum(), v);
314 inline double deltaEta(
const Particle& p,
double eta) {
315 return deltaEta(p.momentum(),
eta);
318 inline double deltaEta(
const FourMomentum& v,
const Particle& p) {
319 return deltaEta(v, p.momentum());
322 inline double deltaEta(
const FourVector& v,
const Particle& p) {
323 return deltaEta(v, p.momentum());
326 inline double deltaEta(
const Vector3& v,
const Particle& p) {
327 return deltaEta(v, p.momentum());
330 inline double deltaEta(
double eta,
const Particle& p) {
331 return deltaEta(eta, p.momentum());
int PdgId
Typedef for a PDG ID code.
Definition: Particle.fhh:29
bool cmpParticleByAscP(const Particle &a, const Particle &b)
Sort by ascending momentum, .
Definition: Particle.hh:166
bool cmpParticleByEt(const Particle &a, const Particle &b)
Sort by descending transverse energy, .
Definition: Particle.hh:170
const FourMomentum & momentum() const
The momentum of this Particle.
Definition: Particle.hh:68
bool cmpParticleByAscAbsRapidity(const Particle &a, const Particle &b)
Sort by ascending absolute rapidity, .
Definition: Particle.hh:214
double E() const
Get energy (time component of momentum).
Definition: Vector4.hh:355
Particle(const GenParticle &gp)
Constructor from a HepMC GenParticle.
Definition: Particle.hh:33
bool cmpParticleByAscRapidity(const Particle &a, const Particle &b)
Sort by ascending rapidity, .
Definition: Particle.hh:206
double mod() const
Calculate the modulus of a vector. .
Definition: VectorN.hh:87
Vector3 vector3() const
Get the spatial part of the 4-vector as a 3-vector.
Definition: Vector4.hh:139
double phi(const Vector3 &v, const PhiMapping mapping=ZERO_2PI)
Synonym for azimuthalAngle.
Definition: Vector3.hh:281
Representation of particles from a HepMC::GenEvent.
Definition: Particle.hh:16
PdgId pdgId() const
The PDG ID code for this Particle.
Definition: Particle.hh:56
double mass() const
The mass of this Particle.
Definition: Particle.hh:78
bool cmpParticleByAscPt(const Particle &a, const Particle &b)
Sort by ascending transverse momentum, .
Definition: Particle.hh:158
double pT() const
Calculate the transverse momentum .
Definition: Vector4.hh:410
bool hasAncestor(PdgId pdg_id) const
Definition: Particle.cc:8
bool cmpParticleByAscE(const Particle &a, const Particle &b)
Sort by ascending energy, .
Definition: Particle.hh:182
std::pair< Particle, Particle > ParticlePair
Typedef for a pair of Particle objects.
Definition: Particle.fhh:20
bool fromDecay() const
Determine whether the particle is from a hadron or tau decay.
Definition: Particle.cc:19
bool cmpParticleByDescRapidity(const Particle &a, const Particle &b)
Sort by descending rapidity, .
Definition: Particle.hh:202
double mass() const
Get the mass (the Lorentz self-invariant).
Definition: Vector4.hh:384
bool cmpParticleByDescAbsRapidity(const Particle &a, const Particle &b)
Sort by descending absolute rapidity, .
Definition: Particle.hh:210
double energy() const
The energy of this Particle.
Definition: Particle.hh:73
bool cmpParticleByAscPseudorapidity(const Particle &a, const Particle &b)
Sort by ascending pseudorapidity, .
Definition: Particle.hh:190
const std::string & toParticleName(PdgId p)
Print a PdgId as a named string.
Definition: ParticleName.hh:142
Base class for particle-like things like Particle and Jet.
Definition: ParticleBase.hh:11
double pseudorapidity() const
Pseudorapidity (defined purely by the 3-vector components)
Definition: Vector4.hh:129
bool cmpParticleByAscEt(const Particle &a, const Particle &b)
Sort by ascending transverse energy, .
Definition: Particle.hh:174
Particle()
Definition: Particle.hh:21
RapScheme
Enum for rapidity variable to be used in calculating , applying rapidity cuts, etc.
Definition: MathHeader.hh:60
const GenParticle & genParticle() const
Get a const reference to the original GenParticle.
Definition: Particle.hh:43
bool cmpParticleByDescAbsPseudorapidity(const Particle &a, const Particle &b)
Sort by descending absolute pseudorapidity, .
Definition: Particle.hh:194
bool cmpParticleByE(const Particle &a, const Particle &b)
Sort by descending energy, .
Definition: Particle.hh:178
bool hasGenParticle() const
Check if the particle corresponds to a GenParticle.
Definition: Particle.hh:50
double Et() const
Calculate the transverse energy .
Definition: Vector4.hh:420
Particle & setMomentum(const FourMomentum &momentum)
Set the momentum of this Particle.
Definition: Particle.hh:62
std::string toString(const AnalysisInfo &ai)
String representation.
Definition: AnalysisInfo.cc:234
double eta(const Vector3 &v)
Synonym for pseudorapidity.
Definition: Vector3.hh:299
bool cmpParticleByP(const Particle &a, const Particle &b)
Sort by descending momentum, .
Definition: Particle.hh:162
bool cmpParticleByPt(const Particle &a, const Particle &b)
Sort by descending transverse momentum, .
Definition: Particle.hh:154
bool cmpParticleByDescPseudorapidity(const Particle &a, const Particle &b)
Sort by descending pseudorapidity, .
Definition: Particle.hh:186
std::ostream & operator<<(std::ostream &os, const AnalysisInfo &ai)
Stream an AnalysisInfo as a text description.
Definition: AnalysisInfo.hh:239
double rapidity() const
Calculate the rapidity.
Definition: Vector4.hh:400
Specialized version of the FourVector with momentum/energy functionality.
Definition: Vector4.hh:324
Particle(PdgId pid, const FourMomentum &mom)
Constructor without GenParticle.
Definition: Particle.hh:27
bool cmpParticleByAscAbsPseudorapidity(const Particle &a, const Particle &b)
Sort by ascending absolute pseudorapidity, .
Definition: Particle.hh:198