Belle II Software development
Particle Class Reference

Class to store reconstructed particles. More...

#include <Particle.h>

Inheritance diagram for Particle:
RelationsInterface< BASE >

Public Types

enum  EParticleSourceObject {
  c_Undefined = 0 ,
  c_Track = 1 ,
  c_ECLCluster = 2 ,
  c_KLMCluster = 3 ,
  c_V0 = 4 ,
  c_MCParticle = 5 ,
  c_Composite = 6 ,
  c_NoMDSTSource = 7
}
 particle source enumerators More...
 
enum  EFlavorType {
  c_Unflavored = 0 ,
  c_Flavored = 1
}
 describes flavor type, see getFlavorType(). More...
 
enum  {
  c_DimPosition = 3 ,
  c_DimMomentum = 4 ,
  c_DimMatrix = 7 ,
  c_SizeMatrix = c_DimMatrix * (c_DimMatrix + 1) / 2
}
 error matrix dimensions and size of 1D representation More...
 
enum  {
  c_Px ,
  c_Py ,
  c_Pz ,
  c_E ,
  c_X ,
  c_Y ,
  c_Z
}
 enumerator used for error matrix handling, shows also in which order the rows (columns) are defined More...
 
enum  PropertyFlags {
  c_Ordinary = 0 ,
  c_IsUnspecified = 1 ,
  c_IsIgnoreRadiatedPhotons = 2 ,
  c_IsIgnoreIntermediate = 4 ,
  c_IsIgnoreMassive = 8 ,
  c_IsIgnoreNeutrino = 16 ,
  c_IsIgnoreGamma = 32 ,
  c_IsIgnoreBrems = 64 ,
  c_IsIgnoreMisID = 128 ,
  c_IsIgnoreDecayInFlight = 256
}
 Flags that describe the particle property, which are used in the MC matching. More...
 

Public Member Functions

 Particle ()
 Default constructor.
 
 Particle (const ROOT::Math::PxPyPzEVector &momentum, const int pdgCode)
 Constructor from a Lorentz vector and PDG code.
 
 Particle (const ROOT::Math::PxPyPzEVector &momentum, const int pdgCode, EFlavorType flavorType, const EParticleSourceObject particleType, const unsigned mdstIndex)
 Constructor for final state particles.
 
 Particle (const ROOT::Math::PxPyPzEVector &momentum, const int pdgCode, EFlavorType flavorType, const std::vector< int > &daughterIndices, TClonesArray *arrayPointer=nullptr)
 Constructor for composite particles.
 
 Particle (const ROOT::Math::PxPyPzEVector &momentum, const int pdgCode, EFlavorType flavorType, const std::vector< int > &daughterIndices, int properties, TClonesArray *arrayPointer=nullptr)
 Constructor for composite particles.
 
 Particle (const ROOT::Math::PxPyPzEVector &momentum, const int pdgCode, EFlavorType flavorType, const std::vector< int > &daughterIndices, int properties, const std::vector< int > &daughterProperties, TClonesArray *arrayPointer=nullptr)
 Constructor for composite particles.
 
 Particle (const Track *track, const Const::ChargedStable &chargedStable)
 Constructor from a reconstructed track (mdst object Track);.
 
 Particle (int trackArrayIndex, const TrackFitResult *trackFit, const Const::ChargedStable &chargedStable)
 Constructor from a reconstructed Track given as TrackFitResult.
 
 Particle (const ECLCluster *eclCluster, const Const::ParticleType &type=Const::photon)
 Constructor of a photon from a reconstructed ECL cluster that is not matched to any charged track.
 
 Particle (const KLMCluster *klmCluster, const int pdgCode=Const::Klong.getPDGCode())
 Constructor from a reconstructed KLM cluster.
 
 Particle (const MCParticle *MCparticle)
 Constructor from MC particle (mdst object MCParticle)
 
 ~Particle ()
 Destructor.
 
void setPDGCode (const int pdg)
 Sets PDG code.
 
void set4Vector (const ROOT::Math::PxPyPzEVector &p4)
 Sets Lorentz vector.
 
void set4VectorDividingByMomentumScaling (const ROOT::Math::PxPyPzEVector &p4)
 Sets Lorentz vector dividing by the momentum scaling factor.
 
void setVertex (const ROOT::Math::XYZVector &vertex)
 Sets position (decay vertex)
 
void setEnergyLossCorrection (double energyLossCorrection)
 Sets Energy loss correction.
 
double getMomentumLossCorrectionFactor () const
 Returns effect of energy correction on the particle momentum.
 
void setMomentumScalingFactor (double momentumScalingFactor)
 Sets momentum scaling.
 
void setMomentumSmearingFactor (double momentumSmearingFactor)
 Sets momentum smearing.
 
void setJacobiMatrix (const TMatrixF &jacobiMatrix)
 Sets 4x6 jacobi matrix.
 
void setMomentumVertexErrorMatrix (const TMatrixFSym &errMatrix)
 Sets 7x7 error matrix.
 
void setPValue (double pValue)
 Sets chi^2 probability of fit.
 
void setProperty (const int properties)
 sets m_properties
 
void updateMomentum (const ROOT::Math::PxPyPzEVector &p4, const ROOT::Math::XYZVector &vertex, const TMatrixFSym &errMatrix, double pValue)
 Sets Lorentz vector, position, 7x7 error matrix and p-value.
 
void updateMass (const int pdgCode)
 Updates particle mass with the mass of the particle corresponding to the given PDG.
 
void appendDaughter (const Particle *daughter, const bool updateType=true, const int daughterProperty=c_Ordinary)
 Appends index of daughter to daughters index array.
 
void appendDaughter (int particleIndex, const bool updateType=true)
 Appends index of daughter to daughters index array.
 
void removeDaughter (const Particle *daughter, const bool updateType=true)
 Removes index of daughter from daughters index array.
 
bool replaceDaughter (const Particle *oldDaughter, Particle *newDaughter)
 Replace index of given daughter with new daughter, return true if a replacement is made.
 
bool replaceDaughterRecursively (const Particle *oldDaughter, Particle *newDaughter)
 Apply replaceDaughter to all Particles in the decay tree by looping recursively through it, return true if a replacement is made.
 
int getPDGCode (void) const
 Returns PDG code.
 
double getCharge (void) const
 Returns particle charge.
 
EFlavorType getFlavorType () const
 Returns flavor type of the decay (for FS particles: flavor type of particle)
 
EParticleSourceObject getParticleSource () const
 Returns particle source as defined with enum EParticleSourceObject.
 
unsigned getMdstArrayIndex (void) const
 Returns 0-based index of MDST store array object (0 for composite particles)
 
int getProperty () const
 Returns particle property as a bit pattern The values are defined in the PropertyFlags enum and described in detail there.
 
double getMass () const
 Returns invariant mass (= nominal for FS particles)
 
double getPDGMass (void) const
 Returns uncertainty on the invariant mass (requires valid momentum error matrix)
 
double getPDGLifetime () const
 Returns particle nominal lifetime.
 
double getEnergy () const
 Returns total energy.
 
ROOT::Math::PxPyPzEVector get4Vector () const
 Returns Lorentz vector.
 
ROOT::Math::XYZVector getMomentum () const
 Returns momentum vector.
 
double getMomentumMagnitude () const
 Returns momentum magnitude.
 
double getP () const
 Returns momentum magnitude (same as getMomentumMagnitude but with shorter name)
 
double getPx () const
 Returns x component of momentum.
 
double getPy () const
 Returns y component of momentum.
 
double getPz () const
 Returns z component of momentum.
 
double getEffectiveMomentumScale () const
 Returns effective momentum scale which is the product of the momentum scaling and smearing factors.
 
double getEnergyLossCorrection () const
 Returns Energy Loss Correction.
 
ROOT::Math::XYZVector getVertex () const
 Returns vertex position (POCA for charged, IP for neutral FS particles)
 
double getX () const
 Returns x component of vertex position.
 
double getY () const
 Returns y component of vertex position.
 
double getZ () const
 Returns z component of vertex position.
 
double getPValue () const
 Returns chi^2 probability of fit if done or -1.
 
TMatrixFSym getMomentumVertexErrorMatrix () const
 Returns 7x7 error matrix.
 
TMatrixFSym getMomentumErrorMatrix () const
 Returns the 4x4 momentum error matrix.
 
TMatrixFSym getVertexErrorMatrix () const
 Returns the 3x3 position error sub-matrix.
 
double getCosHelicity (const Particle *mother=nullptr) const
 Returns cosine of the helicity angle The helicity angle is defined in the rest frame of the particle as the angle between the negative momentum of the mother and.
 
double getCosHelicityDaughter (unsigned iDaughter, unsigned iGrandDaughter=0) const
 Returns cosine of the helicity angle of the given daughter defined by given grand daughter.
 
double getAcoplanarity () const
 Returns acoplanarity angle defined as the angle between the decay planes of the grand daughters in the particle's rest frame This assumes that the particle and its daughters have two daughters each.
 
int getMdstSource () const
 Returns unique identifier of final state particle (needed in particle combiner)
 
unsigned getNDaughters (void) const
 Returns number of daughter particles.
 
const std::vector< int > & getDaughterIndices () const
 Returns a vector of store array indices of daughter particles.
 
const std::vector< int > & getDaughterProperties () const
 Returns a vector of properties of daughter particles.
 
const ParticlegetDaughter (unsigned i) const
 Returns a pointer to the i-th daughter particle.
 
bool forEachDaughter (const std::function< bool(const Particle *)> &function, bool recursive=true, bool includeSelf=true) const
 Apply a function to all daughters of this particle.
 
std::vector< Belle2::Particle * > getDaughters () const
 Returns a vector of pointers to daughter particles.
 
std::vector< const Belle2::Particle * > getFinalStateDaughters () const
 Returns a vector of pointers to Final State daughter particles.
 
std::vector< const Belle2::Particle * > getAllDaughters () const
 Returns a vector of pointers to all generations' daughter particles.
 
std::vector< int > getMdstArrayIndices (EParticleSourceObject type) const
 Returns a vector of StoreArray indices of given MDST dataobjects.
 
bool overlapsWith (const Particle *oParticle) const
 Returns true if final state ancestors of oParticle overlap.
 
bool isCopyOf (const Particle *oParticle, bool doDetailedComparison=false) const
 Returns true if this Particle and oParticle are copies of each other.
 
const TrackgetTrack () const
 Returns the pointer to the Track object that was used to create this Particle (ParticleType == c_Track).
 
const TrackFitResultgetTrackFitResult () const
 Returns the pointer to the TrackFitResult that was used to create this Particle (ParticleType == c_Track).
 
const V0getV0 () const
 Returns the pointer to the V0 object that was used to create this Particle (if ParticleType == c_V0).
 
const PIDLikelihoodgetPIDLikelihood () const
 Returns the pointer to the PIDLikelihood object that is related to the Track, which was used to create this Particle (ParticleType == c_Track).
 
const ECLClustergetECLCluster () const
 Returns the pointer to the ECLCluster object that was used to create this Particle (if ParticleType == c_ECLCluster).
 
double getECLClusterEnergy () const
 Returns the energy of the ECLCluster for the particle.
 
const KLMClustergetKLMCluster () const
 Returns the pointer to the KLMCluster object that was used to create this Particle (ParticleType == c_KLMCluster).
 
const MCParticlegetMCParticle () const
 Returns the pointer to the MCParticle object that was used to create this Particle (ParticleType == c_MCParticle).
 
std::string getName () const override
 Return name of this particle.
 
std::string getInfoHTML () const override
 Return a short summary of this object's contents in HTML format.
 
void print () const
 Prints the contents of a Particle object to standard output.
 
std::vector< std::string > getExtraInfoNames () const
 get a list of the extra info names
 
void removeExtraInfo ()
 Remove all stored extra info fields.
 
double getExtraInfo (const std::string &name) const
 Return given value if set.
 
bool hasExtraInfo (const std::string &name) const
 Return whether the extra info with the given name is set.
 
int getExtraInfoMap () const
 Return the id of the associated ParticleExtraInfoMap or -1 if no map is set.
 
unsigned int getExtraInfoSize () const
 Return the size of the extra info array.
 
void writeExtraInfo (const std::string &name, const double value)
 Sets the user defined extraInfo.
 
void setExtraInfo (const std::string &name, double value)
 Sets the user-defined data of given name to the given value.
 
void addExtraInfo (const std::string &name, double value)
 Sets the user-defined data of given name to the given value.
 
TClonesArray * getArrayPointer () const
 Returns the pointer to the store array which holds the daughter particles.
 
int getPDGCodeUsedForFit () const
 Return the always positive PDG code which was used for the track fit (if there was a track fit) of this particle.
 
bool wasExactFitHypothesisUsed () const
 Returns true if the type represented by this Particle object was used use as a mass hypothesis during the track of this Particle's parameters.
 
bool isMostLikely () const
 Returns true if the (track-based) particle is created with its most likely mass hypothesis based on PID likelihood.
 
std::pair< Const::ChargedStable, const TrackFitResult * > getMostLikelyTrackFitResult () const
 For a (track-based) particle, returns the charged stable mass hypothesis associated to the most probable TrackFitResult, and the TrackFitResult itself.
 
bool isMostLikelyTrackFitResult () const
 Returns true if the (track-based) particle is created with its most likely mass hypothesis based on TrackFitResult.
 
ECLCluster::EHypothesisBit getECLClusterEHypothesisBit () const
 Returns the ECLCluster EHypothesisBit for this Particle.
 
const ParticlegetParticleFromGeneralizedIndexString (const std::string &generalizedIndex) const
 Explores the decay tree of the particle and returns the (grand^n)daughter identified by a generalized index.
 
void updateJacobiMatrix ()
 Propagate the photon energy scaling to jacobian elements that were calculated using energy.
 
void fillFSPDaughters (std::vector< const Belle2::Particle * > &fspDaughters) const
 Fill final state particle daughters into a vector.
 
void fillAllDaughters (std::vector< const Belle2::Particle * > &allDaughters) const
 Fill all generations' daughters into a vector.
 
void addRelationTo (const RelationsInterface< BASE > *object, float weight=1.0, const std::string &namedRelation="") const
 Add a relation from this object to another object (with caching).
 
void addRelationTo (const TObject *object, float weight=1.0, const std::string &namedRelation="") const
 Add a relation from this object to another object (no caching, can be quite slow).
 
void copyRelations (const RelationsInterface< BASE > *sourceObj)
 Copies all relations of sourceObj (pointing from or to sourceObj) to this object (including weights).
 
template<class TO >
RelationVector< TO > getRelationsTo (const std::string &name="", const std::string &namedRelation="") const
 Get the relations that point from this object to another store array.
 
template<class FROM >
RelationVector< FROM > getRelationsFrom (const std::string &name="", const std::string &namedRelation="") const
 Get the relations that point from another store array to this object.
 
template<class T >
RelationVector< T > getRelationsWith (const std::string &name="", const std::string &namedRelation="") const
 Get the relations between this object and another store array.
 
template<class TO >
TO * getRelatedTo (const std::string &name="", const std::string &namedRelation="") const
 Get the object to which this object has a relation.
 
template<class FROM >
FROM * getRelatedFrom (const std::string &name="", const std::string &namedRelation="") const
 Get the object from which this object has a relation.
 
template<class T >
T * getRelated (const std::string &name="", const std::string &namedRelation="") const
 Get the object to or from which this object has a relation.
 
template<class TO >
std::pair< TO *, float > getRelatedToWithWeight (const std::string &name="", const std::string &namedRelation="") const
 Get first related object & weight of relation pointing to an array.
 
template<class FROM >
std::pair< FROM *, float > getRelatedFromWithWeight (const std::string &name="", const std::string &namedRelation="") const
 Get first related object & weight of relation pointing from an array.
 
template<class T >
std::pair< T *, float > getRelatedWithWeight (const std::string &name="", const std::string &namedRelation="") const
 Get first related object & weight of relation pointing from/to an array.
 
std::string getInfo () const
 Return a short summary of this object's contents in raw text format.
 
std::string getArrayName () const
 Get name of array this object is stored in, or "" if not found.
 
int getArrayIndex () const
 Returns this object's array index (in StoreArray), or -1 if not found.
 

Private Member Functions

void setMomentumPositionErrorMatrix (const TrackFitResult *trackFit)
 Sets the momentum, position and error matrix for this particle (created from charged Track)
 
void resetErrorMatrix ()
 Resets 7x7 error matrix All elements are set to 0.0.
 
void resetJacobiMatrix ()
 Resets 4x6 error matrix All elements are set to 0.0.
 
void storeErrorMatrix (const TMatrixFSym &errMatrix)
 Stores 7x7 error matrix into private member m_errMatrix.
 
void storeJacobiMatrix (const TMatrixF &jacobiMatrix)
 Stores 4x6 Jacobi matrix into private member m_jacobiMatrix.
 
void fillDecayChain (std::vector< int > &decayChain) const
 Fill vector with (PDGCode, MdstSource) pairs for the entire decay chain.
 
void setFlavorType ()
 sets m_flavorType using m_pdgCode
 
void setMdstArrayIndex (const int arrayIndex)
 set mdst array index
 
int generatePDGCodeFromCharge (const int chargedSign, const Const::ChargedStable &chargedStable)
 Generate the PDG code with correct sign, using the charge.
 
 ClassDefOverride (Particle, 17)
 Class to store reconstructed particles.
 
 ClassDef (RelationsInterface, 0)
 defines interface for accessing relations of objects in StoreArray.
 

Private Attributes

int m_pdgCode
 PDG code.
 
int m_pdgCodeUsedForFit = 0
 PDG code used for the track fit.
 
double m_mass
 particle (invariant) mass
 
double m_px
 momentum component x
 
double m_py
 momentum component y
 
double m_pz
 momentum component z
 
double m_momentumScale = 1.0
 effective momentum scale factor
 
double m_momentumScalingFactor = 1.0
 momentum scaling factor
 
double m_momentumSmearingFactor = 1.0
 momentum smearing factor
 
double m_energyLossCorrection = 0.0
 energy loss correction.
 
double m_x
 position component x
 
double m_y
 position component y
 
double m_z
 position component z
 
double m_errMatrix [c_SizeMatrix] = {}
 error matrix (1D representation)
 
double m_jacobiMatrix [c_SizeMatrix] = {}
 error matrix (1D representation)
 
double m_pValue
 chi^2 probability of the fit.
 
std::vector< int > m_daughterIndices
 daughter particle indices
 
EFlavorType m_flavorType
 flavor type.
 
EParticleSourceObject m_particleSource
 (mdst) source of particle
 
unsigned m_mdstIndex
 0-based index of MDST store array object
 
int m_properties
 particle property
 
std::vector< int > m_daughterProperties
 daughter particle properties
 
int m_identifier = -1
 Identifier that can be used to identify whether the particle is unique or is a copy or representation of another.
 
std::vector< double > m_extraInfo
 Stores associated user defined values.
 
TClonesArray * m_arrayPointer
 Internal pointer to DataStore array containing the daughters of this particle.
 
DataStore::StoreEntrym_cacheDataStoreEntry
 Cache of the data store entry to which this object belongs.
 
int m_cacheArrayIndex
 Cache of the index in the TClonesArray to which this object belongs.
 

Friends

class ParticleSubset
 

Detailed Description

Class to store reconstructed particles.

This class is a common representation of all particle types, e.g.:

  • final state particles (FS particles):
    • charged kaons/pions/electrons/muons/protons reconstructed as Track
    • photons reconstructed as ECLCluster (without associated Track)
    • long lived neutral kaons reconstructed in KLM (as KLMCluster without associated Track)
  • composite particles:
    • pre-reconstructed V0 particles: Kshort, Lambda baryon, converted photon, ...
    • reconstructed in decays (via combinations)

Private members are limited to those which completely define the particle and that are common to all particle types. These are:

  • particle mass
    • nominal for FS particles
    • invariant for composite particles
  • momentum vector (px, py, pz)
  • position (x, y, z)
    • POCA for charged FS particles
    • IP for photons, Klong and pi0
    • decay vertex for composite particles
  • 7x7 error matrix (order is: px, py, pz, E, x, y, z)
  • chi^2 probability of the fit (track fit, pi0 mass-constr. fit or vertex fit)
  • PDG code
  • vector of StoreArray<Particle> indices of daughter particles

Additional private members are needed in order to make composite particles (via combinations):

  • mdst index of an object from which the FS particle is created
  • source of the object from which the particle is created (see EParticleSourceObject)
  • flavor type (unflavored/flavored) of a decay or flavor type of FS particle

Finally, it is possible to store user-defined floating-point values using addExtraInfo() and getExtraInfo(), identified by a string key. These are general purpose functions for collecting additional information about reconstructed particles that can not be included into any existing data object. Different Particle objects can store different sets of values, and there are no hard-coded limitations on the number of values stored.

Definition at line 75 of file Particle.h.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum

error matrix dimensions and size of 1D representation

Definition at line 102 of file Particle.h.

102 {c_DimPosition = 3, c_DimMomentum = 4, c_DimMatrix = 7,
103 c_SizeMatrix = c_DimMatrix * (c_DimMatrix + 1) / 2
104 };

◆ anonymous enum

anonymous enum

enumerator used for error matrix handling, shows also in which order the rows (columns) are defined

Definition at line 110 of file Particle.h.

110{c_Px, c_Py, c_Pz, c_E, c_X, c_Y, c_Z};

◆ EFlavorType

describes flavor type, see getFlavorType().

Enumerator
c_Unflavored 

Is its own antiparticle or we don't know whether it is a particle/antiparticle.

c_Flavored 

Is either particle or antiparticle.

Definition at line 94 of file Particle.h.

94 {
95 c_Unflavored = 0,
96 c_Flavored = 1,
97 };
@ c_Unflavored
Is its own antiparticle or we don't know whether it is a particle/antiparticle.
Definition: Particle.h:95
@ c_Flavored
Is either particle or antiparticle.
Definition: Particle.h:96

◆ EParticleSourceObject

particle source enumerators

Definition at line 82 of file Particle.h.

82 {
83 c_Undefined = 0,
84 c_Track = 1,
85 c_ECLCluster = 2,
86 c_KLMCluster = 3,
87 c_V0 = 4,
88 c_MCParticle = 5,
89 c_Composite = 6,
90 c_NoMDSTSource = 7
91 };

◆ PropertyFlags

Flags that describe the particle property, which are used in the MC matching.

Enumerator
c_IsUnspecified 

Ordinary particles.

Is the particle unspecified by marking @ ?

c_IsIgnoreRadiatedPhotons 

Is the particle MC matched with the ignore radiated photon flag set?

c_IsIgnoreIntermediate 

Is the particle MC matched with the ignore intermediate resonances flag set?

c_IsIgnoreMassive 

Is the particle MC matched with the ignore missing massive particle flag set?

c_IsIgnoreNeutrino 

Is the particle MC matched with the ignore missing neutrino flag set?

c_IsIgnoreGamma 

Is the particle MC matched with the ignore missing gamma flag set?

c_IsIgnoreBrems 

Is the particle MC matched with the ignore added Brems gamma flag set?

c_IsIgnoreMisID 

Is the particle MC matched with the ignore MisID flag set?

c_IsIgnoreDecayInFlight 

Is the particle MC matched with the ignore DecayInFlight flag set?

Definition at line 116 of file Particle.h.

116 {
117 c_Ordinary = 0,
118 c_IsUnspecified = 1,
122 c_IsIgnoreNeutrino = 16,
123 c_IsIgnoreGamma = 32,
124 c_IsIgnoreBrems = 64,
125 c_IsIgnoreMisID = 128,
127 };
@ c_IsIgnoreNeutrino
Is the particle MC matched with the ignore missing neutrino flag set?
Definition: Particle.h:122
@ c_IsIgnoreRadiatedPhotons
Is the particle MC matched with the ignore radiated photon flag set?
Definition: Particle.h:119
@ c_IsIgnoreGamma
Is the particle MC matched with the ignore missing gamma flag set?
Definition: Particle.h:123
@ c_IsUnspecified
Ordinary particles.
Definition: Particle.h:118
@ c_IsIgnoreBrems
Is the particle MC matched with the ignore added Brems gamma flag set?
Definition: Particle.h:124
@ c_IsIgnoreDecayInFlight
Is the particle MC matched with the ignore DecayInFlight flag set?
Definition: Particle.h:126
@ c_IsIgnoreMisID
Is the particle MC matched with the ignore MisID flag set?
Definition: Particle.h:125
@ c_IsIgnoreIntermediate
Is the particle MC matched with the ignore intermediate resonances flag set?
Definition: Particle.h:120
@ c_IsIgnoreMassive
Is the particle MC matched with the ignore missing massive particle flag set?
Definition: Particle.h:121

Constructor & Destructor Documentation

◆ Particle() [1/11]

Particle ( )

Default constructor.

All private members are set to 0. Particle type is set to c_Undefined.

Definition at line 45 of file Particle.cc.

45 :
46 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
48 m_arrayPointer(nullptr)
49{
51}
double m_pValue
chi^2 probability of the fit.
Definition: Particle.h:1069
double m_pz
momentum component z
Definition: Particle.h:1059
unsigned m_mdstIndex
0-based index of MDST store array object
Definition: Particle.h:1073
double m_py
momentum component y
Definition: Particle.h:1058
double m_x
position component x
Definition: Particle.h:1064
double m_px
momentum component x
Definition: Particle.h:1057
TClonesArray * m_arrayPointer
Internal pointer to DataStore array containing the daughters of this particle.
Definition: Particle.h:1100
EParticleSourceObject m_particleSource
(mdst) source of particle
Definition: Particle.h:1072
void resetErrorMatrix()
Resets 7x7 error matrix All elements are set to 0.0.
Definition: Particle.cc:1079
EFlavorType m_flavorType
flavor type.
Definition: Particle.h:1071
double m_mass
particle (invariant) mass
Definition: Particle.h:1056
double m_z
position component z
Definition: Particle.h:1066
int m_pdgCode
PDG code.
Definition: Particle.h:1054
double m_y
position component y
Definition: Particle.h:1065
int m_properties
particle property
Definition: Particle.h:1074

◆ Particle() [2/11]

Particle ( const ROOT::Math::PxPyPzEVector &  momentum,
const int  pdgCode 
)

Constructor from a Lorentz vector and PDG code.

All other private members are set to their default values (0).

Parameters
momentumLorentz vector
pdgCodePDG code

Definition at line 54 of file Particle.cc.

54 :
55 m_pdgCode(pdgCode), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
57{
59 set4Vector(momentum);
61 // set mass of stable charged particle to its nominal value
64 }
65}
The ParticleType class for identifying different particle types.
Definition: Const.h:408
static const ParticleSet chargedStableSet
set of charged stable particles
Definition: Const.h:618
void setFlavorType()
sets m_flavorType using m_pdgCode
Definition: Particle.cc:1150
void set4Vector(const ROOT::Math::PxPyPzEVector &p4)
Sets Lorentz vector.
Definition: Particle.h:271
void updateMass(const int pdgCode)
Updates particle mass with the mass of the particle corresponding to the given PDG.
Definition: Particle.cc:597

◆ Particle() [3/11]

Particle ( const ROOT::Math::PxPyPzEVector &  momentum,
const int  pdgCode,
EFlavorType  flavorType,
const EParticleSourceObject  particleType,
const unsigned  mdstIndex 
)

Constructor for final state particles.

All other private members are set to their default values (0).

Parameters
momentumLorentz vector
pdgCodePDG code
flavorTypeflavor type
particleTypeparticle source
mdstIndexmdst index

Definition at line 68 of file Particle.cc.

72 :
73 m_pdgCode(pdgCode), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
74 m_pValue(-1), m_flavorType(flavorType), m_particleSource(source), m_properties(0), m_arrayPointer(nullptr)
75{
76 if (flavorType == c_Unflavored and pdgCode < 0)
77 m_pdgCode = -pdgCode;
78
79 setMdstArrayIndex(mdstIndex);
80 set4Vector(momentum);
82 // set mass of stable charged particle to its nominal value
85 }
86}
void setMdstArrayIndex(const int arrayIndex)
set mdst array index
Definition: Particle.cc:350

◆ Particle() [4/11]

Particle ( const ROOT::Math::PxPyPzEVector &  momentum,
const int  pdgCode,
EFlavorType  flavorType,
const std::vector< int > &  daughterIndices,
TClonesArray *  arrayPointer = nullptr 
)

Constructor for composite particles.

All other private members are set to their default values (0).

Parameters
momentumLorentz vector
pdgCodePDG code
flavorTypedecay flavor type
daughterIndicesindices of daughters in StoreArray<Particle>
arrayPointerpointer to store array which stores the daughters, if the particle itself is stored in the same array the pointer can be automatically determined

Definition at line 89 of file Particle.cc.

93 :
94 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
95 m_pValue(-1),
96 m_daughterIndices(daughterIndices),
98 m_properties(0), m_arrayPointer(arrayPointer)
99{
100 m_pdgCode = pdgCode;
101 m_flavorType = flavorType;
102 if (flavorType == c_Unflavored and pdgCode < 0)
103 m_pdgCode = -pdgCode;
104 set4Vector(momentum);
106 // set mass of stable charged particle to its nominal value
109 }
110
111 if (!daughterIndices.empty()) {
112 m_particleSource = c_Composite;
113 if (getArrayPointer() == nullptr) {
114 B2FATAL("Composite Particle (with daughters) was constructed outside StoreArray without specifying daughter array!");
115 }
116 for (unsigned int i = 0; i < m_daughterIndices.size(); i++) {
117 m_daughterProperties.push_back(Particle::PropertyFlags::c_Ordinary);
118 }
119 }
120}
TClonesArray * getArrayPointer() const
Returns the pointer to the store array which holds the daughter particles.
Definition: Particle.h:953
std::vector< int > m_daughterProperties
daughter particle properties
Definition: Particle.h:1075
std::vector< int > m_daughterIndices
daughter particle indices
Definition: Particle.h:1070

◆ Particle() [5/11]

Particle ( const ROOT::Math::PxPyPzEVector &  momentum,
const int  pdgCode,
EFlavorType  flavorType,
const std::vector< int > &  daughterIndices,
int  properties,
TClonesArray *  arrayPointer = nullptr 
)

Constructor for composite particles.

All other private members are set to their default values (0).

Parameters
momentumLorentz vector
pdgCodePDG code
flavorTypedecay flavor type
daughterIndicesindices of daughters in StoreArray<Particle>
propertiesparticle property
arrayPointerpointer to store array which stores the daughters, if the particle itself is stored in the same array the pointer can be automatically determined

Definition at line 122 of file Particle.cc.

127 :
128 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
129 m_pValue(-1),
130 m_daughterIndices(daughterIndices),
132 m_arrayPointer(arrayPointer)
133{
134 m_pdgCode = pdgCode;
135 m_flavorType = flavorType;
136 if (flavorType == c_Unflavored and pdgCode < 0)
137 m_pdgCode = -pdgCode;
138 set4Vector(momentum);
140 // set mass of stable charged particle to its nominal value
143 }
144 m_properties = properties;
145
146 if (!daughterIndices.empty()) {
147 m_particleSource = c_Composite;
148 if (getArrayPointer() == nullptr) {
149 B2FATAL("Composite Particle (with daughters) was constructed outside StoreArray without specifying daughter array!");
150 }
151 for (unsigned int i = 0; i < m_daughterIndices.size(); i++) {
152 m_daughterProperties.push_back(Particle::PropertyFlags::c_Ordinary);
153 }
154 }
155}

◆ Particle() [6/11]

Particle ( const ROOT::Math::PxPyPzEVector &  momentum,
const int  pdgCode,
EFlavorType  flavorType,
const std::vector< int > &  daughterIndices,
int  properties,
const std::vector< int > &  daughterProperties,
TClonesArray *  arrayPointer = nullptr 
)

Constructor for composite particles.

All other private members are set to their default values (0).

Parameters
momentumLorentz vector
pdgCodePDG code
flavorTypedecay flavor type
daughterIndicesindices of daughters in StoreArray<Particle>
propertiesparticle property
daughterPropertiesdaughter particle properties
arrayPointerpointer to store array which stores the daughters, if the particle itself is stored in the same array the pointer can be automatically determined

Definition at line 158 of file Particle.cc.

164 :
165 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
166 m_pValue(-1),
167 m_daughterIndices(daughterIndices),
169 m_daughterProperties(daughterProperties),
170 m_arrayPointer(arrayPointer)
171{
172 m_pdgCode = pdgCode;
173 m_flavorType = flavorType;
174 if (flavorType == c_Unflavored and pdgCode < 0)
175 m_pdgCode = -pdgCode;
176 set4Vector(momentum);
178 // set mass of stable charged particle to its nominal value
181 }
182 m_properties = properties;
183
184 if (!daughterIndices.empty()) {
185 m_particleSource = c_Composite;
186 if (getArrayPointer() == nullptr) {
187 B2FATAL("Composite Particle (with daughters) was constructed outside StoreArray without specifying daughter array!");
188 }
189 }
190}

◆ Particle() [7/11]

Particle ( const Track track,
const Const::ChargedStable chargedStable 
)

Constructor from a reconstructed track (mdst object Track);.

Parameters
trackpointer to Track object
chargedStableType of charged particle

Definition at line 193 of file Particle.cc.

194 :
195 Particle(track ? track->getArrayIndex() : 0, track ? track->getTrackFitResultWithClosestMass(chargedStable) : nullptr,
196 chargedStable)
197{
198}
Particle()
Default constructor.
Definition: Particle.cc:45

◆ Particle() [8/11]

Particle ( int  trackArrayIndex,
const TrackFitResult trackFit,
const Const::ChargedStable chargedStable 
)

Constructor from a reconstructed Track given as TrackFitResult.

To be used to create Particle objects from tracks with full control over the hypothesis (e.g. V0 daughters).

Parameters
trackArrayIndextrack StoreArray index
trackFitpointer to TrackFitResult object
chargedStableType of charged particle

Definition at line 200 of file Particle.cc.

202 :
203 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
205{
206 if (!trackFit) return;
207
208 m_flavorType = c_Flavored; //tracks are charged
209 m_particleSource = c_Track;
210
211 setMdstArrayIndex(trackArrayIndex);
212
214 m_pdgCode = generatePDGCodeFromCharge(trackFit->getChargeSign(), chargedStable);
215
216 // set mass
217 if (TDatabasePDG::Instance()->GetParticle(m_pdgCode) == nullptr)
218 B2FATAL("PDG=" << m_pdgCode << " ***code unknown to TDatabasePDG");
219 m_mass = TDatabasePDG::Instance()->GetParticle(m_pdgCode)->Mass() ;
220
221 // set momentum, position and error matrix
223}
int getPDGCode() const
PDG code.
Definition: Const.h:473
void setMomentumPositionErrorMatrix(const TrackFitResult *trackFit)
Sets the momentum, position and error matrix for this particle (created from charged Track)
Definition: Particle.cc:1001
int m_pdgCodeUsedForFit
PDG code used for the track fit.
Definition: Particle.h:1055
int generatePDGCodeFromCharge(const int chargedSign, const Const::ChargedStable &chargedStable)
Generate the PDG code with correct sign, using the charge.
Definition: Particle.cc:1382
short getChargeSign() const
Return track charge (1 or -1).
Const::ParticleType getParticleType() const
Getter for ParticleType of the mass hypothesis of the track fit.

◆ Particle() [9/11]

Particle ( const ECLCluster eclCluster,
const Const::ParticleType type = Const::photon 
)
explicit

Constructor of a photon from a reconstructed ECL cluster that is not matched to any charged track.

Parameters
eclClusterpointer to ECLCluster object
typethe kind of ParticleType we want (photon by default)

Definition at line 225 of file Particle.cc.

225 :
226 m_pdgCode(type.getPDGCode()), m_mass(type.getMass()), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
228{
229 if (!eclCluster) return;
230
232
233 // returns default vertex from clusterutils (from beam parameters if available)
234 // leave it like that for the moment to make it transparent
235 ClusterUtils C;
236 const XYZVector clustervertex = C.GetIPPosition();
237 setVertex(clustervertex);
238
239 const PxPyPzEVector clustermom = C.Get4MomentumFromCluster(eclCluster, clustervertex, getECLClusterEHypothesisBit());
240 m_px = clustermom.Px();
241 m_py = clustermom.Py();
242 m_pz = clustermom.Pz();
243
244 m_particleSource = c_ECLCluster;
245 setMdstArrayIndex(eclCluster->getArrayIndex());
246
247 // set Chi^2 probability:
248 //TODO: gamma quality can be written here
249 m_pValue = 1;
250
251 // get error matrix.
253
254}
Class to provide momentum-related information from ECLClusters.
Definition: ClusterUtils.h:36
const ROOT::Math::PxPyPzEVector Get4MomentumFromCluster(const ECLCluster *cluster, ECLCluster::EHypothesisBit hypo)
Returns four momentum vector.
Definition: ClusterUtils.cc:25
const ROOT::Math::XYZVector GetIPPosition()
Returns default IP position from beam parameters.
void updateJacobiMatrix()
Propagate the photon energy scaling to jacobian elements that were calculated using energy.
Definition: Particle.cc:257
void setVertex(const ROOT::Math::XYZVector &vertex)
Sets position (decay vertex)
Definition: Particle.h:295
ECLCluster::EHypothesisBit getECLClusterEHypothesisBit() const
Returns the ECLCluster EHypothesisBit for this Particle.
Definition: Particle.h:1001
int getArrayIndex() const
Returns this object's array index (in StoreArray), or -1 if not found.

◆ Particle() [10/11]

Particle ( const KLMCluster klmCluster,
const int  pdgCode = Const::Klong.getPDGCode() 
)
explicit

Constructor from a reconstructed KLM cluster.

Parameters
klmClusterpointer to KLMCluster object
pdgCodePDG code (Klong by default)

Definition at line 295 of file Particle.cc.

295 :
296 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
298{
299 if (!klmCluster) return;
300
301 m_pdgCode = pdgCode;
303
304 set4Vector(klmCluster->getMomentum());
305 setVertex(XYZVector(0, 0, 0)); // so far KLMCluster don't provide reliable / usable position information
306 updateMass(m_pdgCode); // KLMCluster internally use Klong mass, overwrite here to allow neutrons
307
308 m_particleSource = c_KLMCluster;
309 setMdstArrayIndex(klmCluster->getArrayIndex());
310
311 // set Chi^2 probability:
312 //TODO: KL quality can be written here
313 m_pValue = -1;
314
315 // TODO: set error matrix
317 //storeErrorMatrix(klmCluster->getErrorMatrix());
318}
ROOT::Math::PxPyPzEVector getMomentum() const
Get momentum.
Definition: KLMCluster.cc:49

◆ Particle() [11/11]

Particle ( const MCParticle MCparticle)
explicit

Constructor from MC particle (mdst object MCParticle)

Parameters
MCparticlepointer to MCParticle object

Definition at line 321 of file Particle.cc.

321 :
322 m_pdgCode(0), m_mass(0), m_px(0), m_py(0), m_pz(0), m_x(0), m_y(0), m_z(0),
324{
325 if (!mcParticle) return;
326
327 m_pdgCode = mcParticle->getPDG();
328 m_particleSource = c_MCParticle; // TODO: what about daughters if not FS particle?
329
330 setMdstArrayIndex(mcParticle->getArrayIndex());
331
332
334
335 // mass and momentum
336 m_mass = mcParticle->getMass();
337 m_px = mcParticle->getMomentum().X();
338 m_py = mcParticle->getMomentum().Y();
339 m_pz = mcParticle->getMomentum().Z();
340 // production vertex
341 // TODO: good only for FS particles, for composite we must use decay vertex
342 setVertex(mcParticle->getVertex());
343
345}

Member Function Documentation

◆ addExtraInfo()

void addExtraInfo ( const std::string &  name,
double  value 
)

Sets the user-defined data of given name to the given value.

throws std::runtime_error if variable is already set.

Definition at line 1336 of file Particle.cc.

1337{
1338 if (hasExtraInfo(name))
1339 throw std::runtime_error(std::string("addExtraInfo: Value '") + name + "' already set!");
1340
1342 if (!extraInfoMap)
1343 extraInfoMap.create();
1344 if (m_extraInfo.empty()) {
1345 unsigned int mapID = extraInfoMap->getMapForNewVar(name);
1346 m_extraInfo.push_back(mapID);
1347 m_extraInfo.push_back(value);
1348 } else {
1349 unsigned int oldMapID = m_extraInfo[0];
1350 unsigned int insertIndex = m_extraInfo.size();
1351 unsigned int mapID = extraInfoMap->getMapForNewVar(name, oldMapID, insertIndex);
1352
1353 m_extraInfo[0] = mapID; //update map
1354 m_extraInfo.push_back(value); //add value
1355 }
1356}
std::vector< double > m_extraInfo
Stores associated user defined values.
Definition: Particle.h:1091
bool hasExtraInfo(const std::string &name) const
Return whether the extra info with the given name is set.
Definition: Particle.cc:1266
bool create(bool replace=false)
Create a default object in the data store.
Type-safe access to single objects in the data store.
Definition: StoreObjPtr.h:96

◆ addRelationTo() [1/2]

void addRelationTo ( const RelationsInterface< BASE > *  object,
float  weight = 1.0,
const std::string &  namedRelation = "" 
) const
inlineinherited

Add a relation from this object to another object (with caching).

Parameters
objectThe object to which the relation should point.
weightThe weight of the relation.
namedRelationAdditional name for the relation, or "" for the default naming

Definition at line 142 of file RelationsObject.h.

143 {
144 if (object)
146 object, object->m_cacheDataStoreEntry, object->m_cacheArrayIndex, weight, namedRelation);
147 }
void addRelation(const TObject *fromObject, StoreEntry *&fromEntry, int &fromIndex, const TObject *toObject, StoreEntry *&toEntry, int &toIndex, float weight, const std::string &namedRelation)
Add a relation from an object in a store array to another object in a store array.
Definition: DataStore.cc:492
static DataStore & Instance()
Instance of singleton Store.
Definition: DataStore.cc:54
DataStore::StoreEntry * m_cacheDataStoreEntry
Cache of the data store entry to which this object belongs.
int m_cacheArrayIndex
Cache of the index in the TClonesArray to which this object belongs.

◆ addRelationTo() [2/2]

void addRelationTo ( const TObject *  object,
float  weight = 1.0,
const std::string &  namedRelation = "" 
) const
inlineinherited

Add a relation from this object to another object (no caching, can be quite slow).

Parameters
objectThe object to which the relation should point.
weightThe weight of the relation.
namedRelationAdditional name for the relation, or "" for the default naming

Definition at line 155 of file RelationsObject.h.

156 {
157 StoreEntry* toEntry = nullptr;
158 int toIndex = -1;
159 DataStore::Instance().addRelation(this, m_cacheDataStoreEntry, m_cacheArrayIndex, object, toEntry, toIndex, weight, namedRelation);
160 }

◆ appendDaughter() [1/2]

void appendDaughter ( const Particle daughter,
const bool  updateType = true,
const int  daughterProperty = c_Ordinary 
)

Appends index of daughter to daughters index array.

Parameters
daughterpointer to the daughter particle
updateTypebool to set whether particle type should be updated
daughterPropertyproperty of the daughter particle

Definition at line 676 of file Particle.cc.

677{
678 if (updateType) {
679 // is it a composite particle or fsr corrected?
680 m_particleSource = c_Composite;
681 }
682
683 // add daughter index
684 m_daughterIndices.push_back(daughter->getArrayIndex());
685 m_daughterProperties.push_back(daughterProperty);
686}

◆ appendDaughter() [2/2]

void appendDaughter ( int  particleIndex,
const bool  updateType = true 
)
inline

Appends index of daughter to daughters index array.

Parameters
particleIndexindex of daughter in StoreArray<Particle>
updateTypebool to set whether particle type should be updated

Definition at line 416 of file Particle.h.

417 {
418 if (updateType) {
419 // is it a composite particle or fsr corrected?
420 m_particleSource = c_Composite;
421 }
422 m_daughterIndices.push_back(particleIndex);
423 m_daughterProperties.push_back(c_Ordinary);
424 }

◆ copyRelations()

void copyRelations ( const RelationsInterface< BASE > *  sourceObj)
inlineinherited

Copies all relations of sourceObj (pointing from or to sourceObj) to this object (including weights).

Useful if you want to make a complete copy of a StoreArray object to make modifications to it, but retain all information on linked objects.

Note: this only works if sourceObj inherits from the same base (e.g. RelationsObject), and only for related objects that also inherit from the same base.

Definition at line 170 of file RelationsObject.h.

171 {
172 if (!sourceObj)
173 return;
174 auto fromRels = sourceObj->getRelationsFrom<RelationsInterface<BASE>>("ALL");
175 for (unsigned int iRel = 0; iRel < fromRels.size(); iRel++) {
176 fromRels.object(iRel)->addRelationTo(this, fromRels.weight(iRel));
177 }
178
179 auto toRels = sourceObj->getRelationsTo<RelationsInterface<BASE>>("ALL");
180 for (unsigned int iRel = 0; iRel < toRels.size(); iRel++) {
181 this->addRelationTo(toRels.object(iRel), toRels.weight(iRel));
182 }
183 }
void addRelationTo(const RelationsInterface< BASE > *object, float weight=1.0, const std::string &namedRelation="") const
Add a relation from this object to another object (with caching).

◆ fillAllDaughters()

void fillAllDaughters ( std::vector< const Belle2::Particle * > &  allDaughters) const

Fill all generations' daughters into a vector.

Function is called recursively

Parameters
allDaughtersvector of daughter particles

Definition at line 1127 of file Particle.cc.

1128{
1129 // this is FSP
1130 if (getNDaughters() == 0)
1131 return;
1132
1133 // this is not FSP (fill it and go one level down)
1134 for (unsigned i = 0; i < getNDaughters(); i++) {
1135 allDaughters.push_back(getDaughter(i));
1136 getDaughter(i)->fillAllDaughters(allDaughters);
1137 }
1138}
void fillAllDaughters(std::vector< const Belle2::Particle * > &allDaughters) const
Fill all generations' daughters into a vector.
Definition: Particle.cc:1127
unsigned getNDaughters(void) const
Returns number of daughter particles.
Definition: Particle.h:727
const Particle * getDaughter(unsigned i) const
Returns a pointer to the i-th daughter particle.
Definition: Particle.cc:631

◆ fillDecayChain()

void fillDecayChain ( std::vector< int > &  decayChain) const
private

Fill vector with (PDGCode, MdstSource) pairs for the entire decay chain.

Used to determine if two Particles are copies or not.

Function is called recursively

Parameters
decayChainvector of (PDGCode, MdstSource) pairs for each particle in the decay chain of this particle

Definition at line 1140 of file Particle.cc.

1141{
1142 decayChain.push_back(m_pdgCode);
1143 decayChain.push_back(getMdstSource());
1144
1145 for (unsigned i = 0; i < getNDaughters(); i++)
1146 getDaughter(i)->fillDecayChain(decayChain);
1147}
int getMdstSource() const
Returns unique identifier of final state particle (needed in particle combiner)
Definition: Particle.cc:369
void fillDecayChain(std::vector< int > &decayChain) const
Fill vector with (PDGCode, MdstSource) pairs for the entire decay chain.
Definition: Particle.cc:1140

◆ fillFSPDaughters()

void fillFSPDaughters ( std::vector< const Belle2::Particle * > &  fspDaughters) const

Fill final state particle daughters into a vector.

Function is called recursively

Parameters
fspDaughtersvector of daughter particles

Definition at line 1114 of file Particle.cc.

1115{
1116 // this is FSP
1117 if (getNDaughters() == 0) {
1118 fspDaughters.push_back(this);
1119 return;
1120 }
1121
1122 // this is not FSP (go one level down)
1123 for (unsigned i = 0; i < getNDaughters(); i++)
1124 getDaughter(i)->fillFSPDaughters(fspDaughters);
1125}
void fillFSPDaughters(std::vector< const Belle2::Particle * > &fspDaughters) const
Fill final state particle daughters into a vector.
Definition: Particle.cc:1114

◆ forEachDaughter()

bool forEachDaughter ( const std::function< bool(const Particle *)> &  function,
bool  recursive = true,
bool  includeSelf = true 
) const

Apply a function to all daughters of this particle.

Parameters
functionfunction object to run on each daughter. If this function returns true the processing will be stopped immediately.
recursiveif true go through all daughters of daughters as well
includeSelfif true also apply the function to this particle
Returns
true if the function returned true for any of the particles it was applied to

Definition at line 1358 of file Particle.cc.

1360{
1361 std::queue<const Particle*> qq;
1362 // If we include ourselves add only this, otherwise directly all children
1363 if (includeSelf) {
1364 qq.push(this);
1365 } else {
1366 for (size_t i = 0; i < getNDaughters(); ++i) qq.push(getDaughter(i));
1367 }
1368 // Now just repeat until done: take the child, run the functor, remove the
1369 // child, add all children if needed
1370 while (!qq.empty()) {
1371 const Particle* p = qq.front();
1372 if (function(p)) return true;
1373 qq.pop();
1374 // Add children if we go through all children recursively or if we look at
1375 // the current particle: we always want the direct children.
1376 if (recursive || p == this)
1377 for (size_t i = 0; i < p->getNDaughters(); ++i) qq.push(p->getDaughter(i));
1378 }
1379 return false;
1380}
Class to store reconstructed particles.
Definition: Particle.h:75

◆ generatePDGCodeFromCharge()

int generatePDGCodeFromCharge ( const int  chargedSign,
const Const::ChargedStable chargedStable 
)
private

Generate the PDG code with correct sign, using the charge.

Parameters
chargedSigncharge of the particle
chargedStableType of charged particle

Definition at line 1382 of file Particle.cc.

1383{
1384 int absPDGCode = chargedStable.getPDGCode();
1385 int PDGCode = absPDGCode * chargeSign;
1386 // flip sign of PDG code for leptons: their PDG code is positive if the lepton charge is negative and vice versa
1387 if (chargedStable == Const::muon || chargedStable == Const::electron) PDGCode = -PDGCode;
1388 return PDGCode;
1389}
static const ChargedStable muon
muon particle
Definition: Const.h:660
static const ChargedStable electron
electron particle
Definition: Const.h:659

◆ get4Vector()

ROOT::Math::PxPyPzEVector get4Vector ( ) const
inline

Returns Lorentz vector.

Returns
Lorentz vector

Definition at line 547 of file Particle.h.

548 {
549 double correction = getMomentumLossCorrectionFactor();
550 return ROOT::Math::PxPyPzEVector(m_momentumScale * correction * m_px,
551 m_momentumScale * correction * m_py,
552 m_momentumScale * correction * m_pz,
553 getEnergy());
554 }
double m_momentumScale
effective momentum scale factor
Definition: Particle.h:1060
double getEnergy() const
Returns total energy.
Definition: Particle.h:535
double getMomentumLossCorrectionFactor() const
Returns effect of energy correction on the particle momentum.
Definition: Particle.h:315

◆ getAcoplanarity()

double getAcoplanarity ( ) const

Returns acoplanarity angle defined as the angle between the decay planes of the grand daughters in the particle's rest frame This assumes that the particle and its daughters have two daughters each.

Returns
acoplanarity angle

Definition at line 533 of file Particle.cc.

534{
535 // check that we have a decay to two daughters and then two grand daughters each
536 if (getNDaughters() != 2) {
537 B2ERROR("Cannot calculate acoplanarity of particle 'name' because the number of daughters is not 2"
538 << LogVar("name", getName()) << LogVar("# of daughters", getNDaughters()));
539 return std::numeric_limits<double>::quiet_NaN();
540 }
541 const Particle* daughter0 = getDaughter(0);
542 const Particle* daughter1 = getDaughter(1);
543 if ((daughter0->getNDaughters() != 2) || (daughter1->getNDaughters() != 2)) {
544 B2ERROR("Cannot calculate acoplanarity of particle 'name' because the number of grand daughters is not 2"
545 << LogVar("name", getName()) << LogVar("# of grand daughters of first daughter", daughter0->getNDaughters())
546 << LogVar("# of grand daughters of second daughter", daughter1->getNDaughters()));
547 return std::numeric_limits<double>::quiet_NaN();
548 }
549
550 // boost vector to the rest frame of the particle
551 Boost boost(get4Vector().BoostToCM());
552
553 // momenta of the daughters and grand daughters in the particle's rest frame
554 PxPyPzEVector pDaughter0 = boost * daughter0->get4Vector();
555 PxPyPzEVector pGrandDaughter0 = boost * daughter0->getDaughter(0)->get4Vector();
556 PxPyPzEVector pDaughter1 = boost * daughter1->get4Vector();
557 PxPyPzEVector pGrandDaughter1 = boost * daughter1->getDaughter(0)->get4Vector();
558
559 // calculate angle between normal vectors
560 XYZVector normal0 = pDaughter0.Vect().Cross(pGrandDaughter0.Vect());
561 XYZVector normal1 = -pDaughter1.Vect().Cross(pGrandDaughter1.Vect());
562 double result = acos(normal0.Unit().Dot(normal1.Unit()));
563 if (normal0.Cross(normal1).Dot(pDaughter0.Vect()) < 0) result = -result;
564
565 return result;
566}
std::string getName() const override
Return name of this particle.
Definition: Particle.cc:1165
ROOT::Math::PxPyPzEVector get4Vector() const
Returns Lorentz vector.
Definition: Particle.h:547
Class to store variables with their name which were sent to the logging service.

◆ getAllDaughters()

std::vector< const Belle2::Particle * > getAllDaughters ( ) const

Returns a vector of pointers to all generations' daughter particles.

Returns
vector of pointers to all generations' daughter particles

Definition at line 657 of file Particle.cc.

658{
659 std::vector<const Particle*> allDaughters;
660 fillAllDaughters(allDaughters);
661
662 return allDaughters;
663}

◆ getArrayIndex()

int getArrayIndex ( ) const
inlineinherited

Returns this object's array index (in StoreArray), or -1 if not found.

Definition at line 385 of file RelationsObject.h.

386 {
388 return m_cacheArrayIndex;
389 }
bool findStoreEntry(const TObject *object, StoreEntry *&entry, int &index)
Find an object in an array in the data store.
Definition: DataStore.cc:398

◆ getArrayName()

std::string getArrayName ( ) const
inlineinherited

Get name of array this object is stored in, or "" if not found.

Definition at line 377 of file RelationsObject.h.

◆ getArrayPointer()

TClonesArray * getArrayPointer ( ) const
inline

Returns the pointer to the store array which holds the daughter particles.

Warning
TClonesArray is dangerously easy to misuse, please avoid.

Definition at line 953 of file Particle.h.

954 {
955 if (!m_arrayPointer)
957 return m_arrayPointer;
958 }
TClonesArray * getArrayPointer() const
Returns the pointer to the raw DataStore array holding this object (protected since these arrays are ...

◆ getCharge()

double getCharge ( void  ) const

Returns particle charge.

Returns
particle charge in units of elementary charge

Definition at line 622 of file Particle.cc.

623{
624 if (TDatabasePDG::Instance()->GetParticle(m_pdgCode) == nullptr) {
625 B2ERROR("PDG=" << m_pdgCode << " ***code unknown to TDatabasePDG");
626 return 0.0;
627 }
628 return TDatabasePDG::Instance()->GetParticle(m_pdgCode)->Charge() / 3.0;
629}

◆ getCosHelicity()

double getCosHelicity ( const Particle mother = nullptr) const

Returns cosine of the helicity angle The helicity angle is defined in the rest frame of the particle as the angle between the negative momentum of the mother and.

  • the momentum of the first daughter for two body decays
  • the momentum of the photon for pi0 Dalitz decays
  • the direction perpendicular to the daughter momenta for three body decays
    Parameters
    mothermother particle, if not given the center of mass system is taken as mother frame
    Returns
    cosine of the helicity angle

Definition at line 459 of file Particle.cc.

460{
461 // boost vector to the rest frame of the particle
462 Boost boost(get4Vector().BoostToCM());
463
464 // momentum of the mother in the particle's rest frame
465 PxPyPzEVector pMother;
466 if (mother) {
467 pMother = mother->get4Vector();
468 } else {
469 static DBObjPtr<CollisionBoostVector> cmsBoost;
471 pMother.SetE(cmsMass->getMass());
472 pMother = Boost(cmsBoost->getBoost()) * pMother;
473 }
474 pMother = boost * pMother;
475
476 // momentum of the daughter (or normal vector) in the particle's rest frame
477 PxPyPzEVector pDaughter;
478 if (getNDaughters() == 2) { // two body decay
479 pDaughter = boost * getDaughter(0)->get4Vector();
480 } else if (getNDaughters() == 3) {
481 if (getPDGCode() == Const::pi0.getPDGCode()) { // pi0 Dalitz decay
482 for (auto& daughter : getDaughters()) {
483 if (daughter->getPDGCode() == Const::photon.getPDGCode()) {
484 pDaughter = daughter->get4Vector();
485 }
486 }
487 pDaughter = boost * pDaughter;
488 } else { // three body decay
489 PxPyPzEVector pDaughter0 = boost * getDaughter(0)->get4Vector();
490 PxPyPzEVector pDaughter1 = boost * getDaughter(1)->get4Vector();
491
492 XYZVector pDaughterNormal(pDaughter0.Vect().Cross(pDaughter1.Vect()));
493 pDaughter.SetPxPyPzE(pDaughterNormal.X(), pDaughterNormal.Y(), pDaughterNormal.Z(), 0); // energy doesn't matter
494 }
495 }
496
497 double mag2 = pMother.P2() * pDaughter.P2();
498 if (mag2 <= 0) return std::numeric_limits<double>::quiet_NaN();
499 return -pMother.Vect().Dot(pDaughter.Vect()) / sqrt(mag2);
500}
static const ParticleType pi0
neutral pion particle
Definition: Const.h:674
static const ParticleType photon
photon particle
Definition: Const.h:673
Class for accessing objects in the database.
Definition: DBObjPtr.h:21
int getPDGCode(void) const
Returns PDG code.
Definition: Particle.h:454
std::vector< Belle2::Particle * > getDaughters() const
Returns a vector of pointers to daughter particles.
Definition: Particle.cc:637
double sqrt(double a)
sqrt for double
Definition: beamHelpers.h:28

◆ getCosHelicityDaughter()

double getCosHelicityDaughter ( unsigned  iDaughter,
unsigned  iGrandDaughter = 0 
) const

Returns cosine of the helicity angle of the given daughter defined by given grand daughter.

Parameters
iDaughter0-based index of daughter particle
iGrandDaughter0-based index of grand daughter particle
Returns
cosine of the helicity angle

Definition at line 502 of file Particle.cc.

503{
504 // check existence of daughter
505 if (getNDaughters() <= iDaughter) {
506 B2ERROR("No daughter of particle 'name' with index 'iDaughter' for calculation of helicity angle"
507 << LogVar("name", getName()) << LogVar("iDaughter", iDaughter));
508 return std::numeric_limits<double>::quiet_NaN();
509 }
510
511 // boost vector to the rest frame of the daughter particle
512 const Particle* daughter = getDaughter(iDaughter);
513 Boost boost(daughter->get4Vector().BoostToCM());
514
515 // momentum of the this particle in the daughter's rest frame
516 PxPyPzEVector pMother = boost * get4Vector();
517
518 // check existence of grand daughter
519 if (daughter->getNDaughters() <= iGrandDaughter) {
520 B2ERROR("No grand daughter of daughter 'iDaughter' of particle 'name' with index 'iGrandDaughter' for calculation of helicity angle"
521 << LogVar("name", getName()) << LogVar("iDaughter", iDaughter) << LogVar("iGrandDaughter", iGrandDaughter));
522 return std::numeric_limits<double>::quiet_NaN();
523 }
524
525 // momentum of the grand daughter in the daughter's rest frame
526 PxPyPzEVector pGrandDaughter = boost * daughter->getDaughter(iGrandDaughter)->get4Vector();
527
528 double mag2 = pMother.P2() * pGrandDaughter.P2();
529 if (mag2 <= 0) return std::numeric_limits<double>::quiet_NaN();
530 return -pMother.Vect().Dot(pGrandDaughter.Vect()) / sqrt(mag2);
531}

◆ getDaughter()

const Particle * getDaughter ( unsigned  i) const

Returns a pointer to the i-th daughter particle.

Parameters
i0-based index of daughter particle
Returns
Pointer to i-th daughter particles

Definition at line 631 of file Particle.cc.

632{
633 if (i >= getNDaughters()) return nullptr;
634 return static_cast<Particle*>(getArrayPointer()->At(m_daughterIndices[i]));
635}

◆ getDaughterIndices()

const std::vector< int > & getDaughterIndices ( ) const
inline

Returns a vector of store array indices of daughter particles.

Returns
vector of store array indices of daughter particle

Definition at line 736 of file Particle.h.

737 {
738 return m_daughterIndices;
739 }

◆ getDaughterProperties()

const std::vector< int > & getDaughterProperties ( ) const
inline

Returns a vector of properties of daughter particles.

Returns
vector of daughter properties

Definition at line 745 of file Particle.h.

746 {
748 }

◆ getDaughters()

std::vector< Belle2::Particle * > getDaughters ( ) const

Returns a vector of pointers to daughter particles.

Returns
vector of pointers to daughter particles

Definition at line 637 of file Particle.cc.

638{
639 const unsigned int nDaughters = getNDaughters();
640 std::vector<Particle*> daughters(nDaughters);
641
642 const TClonesArray* array = getArrayPointer();
643 for (unsigned i = 0; i < nDaughters; i++)
644 daughters[i] = static_cast<Particle*>(array->At(m_daughterIndices[i]));
645
646 return daughters;
647}

◆ getECLCluster()

const ECLCluster * getECLCluster ( ) const

Returns the pointer to the ECLCluster object that was used to create this Particle (if ParticleType == c_ECLCluster).

Returns the pointer to the most energetic ECLCluster matched to the track (if ParticleType == c_Track). NULL pointer is returned, if the Particle has no relation to an ECLCluster (either the particle is a different type or there was no track match).

Returns
const pointer to the ECLCluster

Definition at line 891 of file Particle.cc.

892{
893 if (m_particleSource == c_ECLCluster) {
894 StoreArray<ECLCluster> eclClusters{};
895 return eclClusters[m_mdstIndex];
896 } else if (m_particleSource == c_Track) {
897 // a track may be matched to several clusters under different hypotheses
898 // take the most energetic of the c_nPhotons hypothesis as "the" cluster
899 StoreArray<Track> tracks{};
900 const ECLCluster* bestTrackMatchedCluster = nullptr;
901 double highestEnergy = -1.0;
902 // loop over all clusters matched to this track
903 for (const ECLCluster& cluster : tracks[m_mdstIndex]->getRelationsTo<ECLCluster>()) {
904 // ignore everything except the nPhotons hypothesis
905 if (!cluster.hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons))
906 continue;
907 // check if we're most energetic thus far
908 if (cluster.getEnergy(ECLCluster::EHypothesisBit::c_nPhotons) > highestEnergy) {
909 highestEnergy = cluster.getEnergy(ECLCluster::EHypothesisBit::c_nPhotons);
910 bestTrackMatchedCluster = &cluster;
911 }
912 }
913 return bestTrackMatchedCluster;
914 } else {
915 return nullptr;
916 }
917}
ECL cluster data.
Definition: ECLCluster.h:27
@ c_nPhotons
CR is split into n photons (N1)
Accessor to arrays stored in the data store.
Definition: StoreArray.h:113

◆ getECLClusterEHypothesisBit()

ECLCluster::EHypothesisBit getECLClusterEHypothesisBit ( ) const
inline

Returns the ECLCluster EHypothesisBit for this Particle.

Definition at line 1001 of file Particle.h.

1002 {
1003 const int pdg = abs(getPDGCode());
1004 if ((pdg == Const::photon.getPDGCode())
1005 or (pdg == Const::electron.getPDGCode())
1006 or (pdg == Const::muon.getPDGCode())
1007 or (pdg == Const::pion.getPDGCode())
1008 or (pdg == Const::kaon.getPDGCode())
1009 or (pdg == Const::proton.getPDGCode())
1010 or (pdg == Const::deuteron.getPDGCode())) {
1012 } else if ((pdg == Const::Klong.getPDGCode())
1013 or (pdg == Const::neutron.getPDGCode())) {
1015 } else {
1017 }
1018 }
static const ParticleType neutron
neutron particle
Definition: Const.h:675
static const ChargedStable pion
charged pion particle
Definition: Const.h:661
static const ParticleType Klong
K^0_L particle.
Definition: Const.h:678
static const ChargedStable proton
proton particle
Definition: Const.h:663
static const ChargedStable kaon
charged kaon particle
Definition: Const.h:662
static const ChargedStable deuteron
deuteron particle
Definition: Const.h:664
@ c_neutralHadron
CR is reconstructed as a neutral hadron (N2)
@ c_none
None as initializer.

◆ getECLClusterEnergy()

double getECLClusterEnergy ( ) const

Returns the energy of the ECLCluster for the particle.

The return value depends on the ECLCluster hypothesis.

Returns
energy of the ECLCluster

Definition at line 919 of file Particle.cc.

920{
921 const ECLCluster* cluster = this->getECLCluster();
922 if (!cluster) return 0;
923 return cluster->getEnergy(this->getECLClusterEHypothesisBit());
924}
const ECLCluster * getECLCluster() const
Returns the pointer to the ECLCluster object that was used to create this Particle (if ParticleType =...
Definition: Particle.cc:891

◆ getEffectiveMomentumScale()

double getEffectiveMomentumScale ( ) const
inline

Returns effective momentum scale which is the product of the momentum scaling and smearing factors.

Returns
momentum scaling factor

Definition at line 614 of file Particle.h.

615 {
616 return m_momentumScale;
617 }

◆ getEnergy()

double getEnergy ( ) const
inline

Returns total energy.

Returns
total energy

Definition at line 535 of file Particle.h.

536 {
541 }
double m_energyLossCorrection
energy loss correction.
Definition: Particle.h:1063

◆ getEnergyLossCorrection()

double getEnergyLossCorrection ( ) const
inline

Returns Energy Loss Correction.

Returns
Energy Loss Correction

Definition at line 623 of file Particle.h.

624 {
626 }

◆ getExtraInfo()

double getExtraInfo ( const std::string &  name) const

Return given value if set.

throws std::runtime_error if variable is not set.

Definition at line 1289 of file Particle.cc.

1290{
1291 if (m_extraInfo.empty())
1292 throw std::runtime_error(std::string("getExtraInfo: Value '") + name + "' not found in Particle!");
1293
1294 //get index for name
1295 const auto mapID = (unsigned int)m_extraInfo[0];
1297 if (!extraInfoMap) {
1298 B2FATAL("ParticleExtraInfoMap not available, but needed for storing extra info in Particle!");
1299 }
1300 unsigned int index = extraInfoMap->getIndex(mapID, name);
1301 if (index == 0 or index >= m_extraInfo.size()) //actually indices start at 1
1302 throw std::runtime_error(std::string("getExtraInfo: Value '") + name + "' not found in Particle!");
1303
1304 return m_extraInfo[index];
1305
1306}

◆ getExtraInfoMap()

int getExtraInfoMap ( ) const
inline

Return the id of the associated ParticleExtraInfoMap or -1 if no map is set.

Definition at line 921 of file Particle.h.

922 {
923 return !m_extraInfo.empty() ? m_extraInfo[0] : -1;
924 }

◆ getExtraInfoNames()

std::vector< std::string > getExtraInfoNames ( ) const

get a list of the extra info names

Definition at line 1175 of file Particle.cc.

1176{
1177 std::vector<std::string> out;
1178 if (!m_extraInfo.empty()) {
1180 if (!extraInfoMap)
1181 B2FATAL("ParticleExtraInfoMap not available, but needed for storing extra info in Particle!");
1182 const ParticleExtraInfoMap::IndexMap& map = extraInfoMap->getMap(m_extraInfo[0]);
1183 for (auto const& ee : map) out.push_back(ee.first);
1184 }
1185 return out;
1186}
std::map< std::string, unsigned int > IndexMap
string -> index map.

◆ getExtraInfoSize()

unsigned int getExtraInfoSize ( ) const
inline

Return the size of the extra info array.

Definition at line 927 of file Particle.h.

928 {
929 return m_extraInfo.size();
930 }

◆ getFinalStateDaughters()

std::vector< const Belle2::Particle * > getFinalStateDaughters ( ) const

Returns a vector of pointers to Final State daughter particles.

Returns
vector of pointers to final state daughter particles

Definition at line 649 of file Particle.cc.

650{
651 std::vector<const Particle*> fspDaughters;
652 fillFSPDaughters(fspDaughters);
653
654 return fspDaughters;
655}

◆ getFlavorType()

EFlavorType getFlavorType ( ) const
inline

Returns flavor type of the decay (for FS particles: flavor type of particle)

Returns
flavor type (0=unflavored, 1=flavored)

Definition at line 469 of file Particle.h.

470 {
471 return m_flavorType;
472 }

◆ getInfo()

std::string getInfo ( ) const
inlineinherited

Return a short summary of this object's contents in raw text format.

Returns the contents of getInfoHTML() while translating line-breaks etc.

Note
: You don't need to implement this function (it's not virtual), getInfoHTML() is enough.

Definition at line 370 of file RelationsObject.h.

371 {
373 }
virtual std::string getInfoHTML() const
Return a short summary of this object's contents in HTML format.
std::string htmlToPlainText(const std::string &html)
See RelationsObject::getInfo()

◆ getInfoHTML()

std::string getInfoHTML ( ) const
overridevirtual

Return a short summary of this object's contents in HTML format.

Reimplemented from RelationsInterface< BASE >.

Definition at line 1188 of file Particle.cc.

1189{
1190 std::stringstream stream;
1191 stream << std::setprecision(4);
1192 stream << "<b>collection</b>=" << getArrayName();
1193 stream << "<br>";
1194 stream << " <b>PDGCode</b>=" << m_pdgCode;
1195 stream << " <b>Charge</b>=" << getCharge();
1196 stream << " <b>PDGMass</b>=" << getPDGMass();
1197 stream << "<br>";
1198 stream << " <b>flavorType</b>=" << m_flavorType;
1199 stream << " <b>particleSource</b>=" << m_particleSource;
1200 stream << " <b>particleTypeUsedForFit</b>=" << m_pdgCodeUsedForFit;
1201 stream << "<br>";
1202
1203 stream << " <b>mdstIndex</b>=" << m_mdstIndex;
1204 stream << " <b>arrayIndex</b>=" << getArrayIndex();
1205 stream << " <b>identifier</b>=" << m_identifier;
1206 stream << " <b>daughterIndices</b>: ";
1207 for (int daughterIndex : m_daughterIndices) {
1208 stream << daughterIndex << ", ";
1209 }
1210 if (m_daughterIndices.empty()) stream << " (none)";
1211 stream << "<br>";
1212
1213 if (!m_daughterIndices.empty()) {
1214 stream << " <b>daughter PDGCodes</b>: ";
1215 for (unsigned i = 0; i < m_daughterIndices.size(); i++) {
1216 const Particle* p = getDaughter(i);
1217 if (p) {stream << p->getPDGCode() << ", ";}
1218 else {stream << "?, ";}
1219 }
1220 stream << "<br>";
1221 }
1222
1223 stream << " <b>mass</b>=" << m_mass;
1224 stream << "<br>";
1225
1226 stream << " <b>momentum</b>=" << HTML::getString(ROOT::Math::XYZVector(getPx(), getPy(), getPz()));
1227 stream << " <b>p</b>=" << getP();
1228 stream << "<br>";
1229
1230 stream << " <b>momentum scaling factor</b>=" << m_momentumScale;
1231 stream << "<br>";
1232
1233 stream << " <b>Energy loss correction</b>=" << m_energyLossCorrection;
1234 stream << "<br>";
1235
1236 stream << " <b>position</b>=" << HTML::getString(ROOT::Math::XYZVector(m_x, m_y, m_z));
1237 stream << "<br>";
1238
1239 stream << " <b>p-value of fit</b> (if done): ";
1240 stream << m_pValue;
1241 stream << "<br>";
1242
1243 stream << " <b>error matrix</b> (px, py, pz, E, x, y ,z):<br>";
1245
1246 stream << " <b>extra info</b>=( ";
1247 if (!m_extraInfo.empty()) {
1249 if (!extraInfoMap) {
1250 B2FATAL("ParticleExtraInfoMap not available, but needed for storing extra info in Particle!");
1251 }
1252 const ParticleExtraInfoMap::IndexMap& map = extraInfoMap->getMap(m_extraInfo[0]);
1253 const unsigned int nVars = m_extraInfo.size();
1254 for (const auto& pair : map) {
1255 if (pair.second < nVars) {
1256 stream << pair.first << "=" << m_extraInfo[pair.second] << " ";
1257 }
1258 }
1259
1260 }
1261 stream << ") " << "<br>";
1262
1263 return stream.str();
1264}
double getPx() const
Returns x component of momentum.
Definition: Particle.h:587
double getPz() const
Returns z component of momentum.
Definition: Particle.h:605
double getPy() const
Returns y component of momentum.
Definition: Particle.h:596
double getPDGMass(void) const
Returns uncertainty on the invariant mass (requires valid momentum error matrix)
Definition: Particle.cc:604
double getCharge(void) const
Returns particle charge.
Definition: Particle.cc:622
TMatrixFSym getMomentumVertexErrorMatrix() const
Returns 7x7 error matrix.
Definition: Particle.cc:420
double getP() const
Returns momentum magnitude (same as getMomentumMagnitude but with shorter name)
Definition: Particle.h:578
int m_identifier
Identifier that can be used to identify whether the particle is unique or is a copy or representation...
Definition: Particle.h:1084
std::string getArrayName() const
Get name of array this object is stored in, or "" if not found.
std::string getString(const TMatrixFBase &matrix, int precision=2, bool color=true)
get HTML table representing a matrix.
Definition: HTML.cc:24

◆ getKLMCluster()

const KLMCluster * getKLMCluster ( ) const

Returns the pointer to the KLMCluster object that was used to create this Particle (ParticleType == c_KLMCluster).

Returns the pointer to the KLMCluster object associated to this Particle if ParticleType == c_Track. NULL pointer is returned, if the Particle has no relation to the KLMCluster.

Returns
const pointer to the KLMCluster

Definition at line 926 of file Particle.cc.

927{
928 if (m_particleSource == c_KLMCluster) {
929 StoreArray<KLMCluster> klmClusters{};
930 return klmClusters[m_mdstIndex];
931 } else if (m_particleSource == c_Track) {
932 // If there is an associated KLMCluster, it's the closest one
933 StoreArray<Track> tracks{};
934 const KLMCluster* klmCluster = tracks[m_mdstIndex]->getRelatedTo<KLMCluster>();
935 return klmCluster;
936 } else {
937 return nullptr;
938 }
939}
KLM cluster data.
Definition: KLMCluster.h:28

◆ getMass()

double getMass ( ) const
inline

Returns invariant mass (= nominal for FS particles)

Returns
invariant mass

Definition at line 507 of file Particle.h.

508 {
509 return m_mass;
510 }

◆ getMCParticle()

const MCParticle * getMCParticle ( ) const

Returns the pointer to the MCParticle object that was used to create this Particle (ParticleType == c_MCParticle).

Returns the best MC match for this particle (if ParticleType == c_Track or c_ECLCluster or c_KLMCluster) NULL pointer is returned, if the Particle was not made from MCParticle or not matched.

Returns
const pointer to the MCParticle

Definition at line 942 of file Particle.cc.

943{
944 if (m_particleSource == c_MCParticle) {
945 StoreArray<MCParticle> mcParticles{};
946 return mcParticles[m_mdstIndex];
947 } else {
948 const MCParticle* related = this->getRelated<MCParticle>();
949 if (related)
950 return related;
951 }
952 return nullptr;
953}
A Class to store the Monte Carlo particle information.
Definition: MCParticle.h:32

◆ getMdstArrayIndex()

unsigned getMdstArrayIndex ( void  ) const
inline

Returns 0-based index of MDST store array object (0 for composite particles)

Returns
index of MDST store array object

Definition at line 487 of file Particle.h.

488 {
489 return m_mdstIndex;
490 }

◆ getMdstArrayIndices()

std::vector< int > getMdstArrayIndices ( EParticleSourceObject  type) const

Returns a vector of StoreArray indices of given MDST dataobjects.

Parameters
typeEParticleSourceObject corresponding to a given MDST dataobject
Returns
vector of StoreArray indices of a given MDST dataobjects

Definition at line 665 of file Particle.cc.

666{
667 std::vector<int> mdstIndices;
668 for (const Particle* fsp : getFinalStateDaughters()) {
669 // is this FSP daughter constructed from given MDST source
670 if (fsp->getParticleSource() == source)
671 mdstIndices.push_back(fsp->getMdstArrayIndex());
672 }
673 return mdstIndices;
674}
std::vector< const Belle2::Particle * > getFinalStateDaughters() const
Returns a vector of pointers to Final State daughter particles.
Definition: Particle.cc:649

◆ getMdstSource()

int getMdstSource ( ) const

Returns unique identifier of final state particle (needed in particle combiner)

Returns
unique identifier of final state particle

Definition at line 369 of file Particle.cc.

370{
371 // Is identifier already set
372 if (m_identifier > -1)
373 return m_identifier + (m_particleSource << 24);
374
375 // Identifier is not set.
376 int identifier = 0;
377 if (m_particleSource == c_ECLCluster) {
378 const ECLCluster* cluster = this->getECLCluster();
379 if (cluster) {
380 identifier = cluster->getUniqueClusterId();
381 } else {
382 B2ERROR("Particle is of type = ECLCluster has identifier not set and no relation to ECLCluster.\n"
383 "This has happen because old microDST is analysed with newer version of software.");
384 }
385 } else {
386 identifier = m_mdstIndex;
387 }
388
389 return identifier + (m_particleSource << 24);
390}

◆ getMomentum()

ROOT::Math::XYZVector getMomentum ( ) const
inline

Returns momentum vector.

Returns
momentum vector

Definition at line 560 of file Particle.h.

561 {
562 return m_momentumScale * getMomentumLossCorrectionFactor() * ROOT::Math::XYZVector(m_px, m_py, m_pz);
563 };

◆ getMomentumErrorMatrix()

TMatrixFSym getMomentumErrorMatrix ( ) const

Returns the 4x4 momentum error matrix.

Returns
The 4x4 momentum error matrix (order: px,py,pz,E)

Definition at line 435 of file Particle.cc.

436{
437 TMatrixFSym mom;
438 const TMatrixFSym& full = getMomentumVertexErrorMatrix();
439
440 // get 4x4 (momentum) submatrix from the full error matrix
441 // momentum related elements are in [0,...,3]x[0,...,3] block
442 full.GetSub(0, 3, mom, "S");
443
444 return mom;
445}

◆ getMomentumLossCorrectionFactor()

double getMomentumLossCorrectionFactor ( ) const
inline

Returns effect of energy correction on the particle momentum.

Returns
momentum change

Definition at line 315 of file Particle.h.

316 {
317 if (m_energyLossCorrection == 0.0) {
318 return 1.0;
319 }
320 double origP = m_momentumScale * sqrt(m_px * m_px + m_py * m_py + m_pz * m_pz);
321 if (origP == 0.0) {
322 return 1.0;
323 }
324 double origE = sqrt(origP * origP + m_mass * m_mass);
325
326 double newP = sqrt((origE - m_energyLossCorrection) * (origE - m_energyLossCorrection) - m_mass * m_mass);
327 return newP / origP;
328 }

◆ getMomentumMagnitude()

double getMomentumMagnitude ( ) const
inline

Returns momentum magnitude.

Returns
momentum magnitude

Definition at line 569 of file Particle.h.

570 {
571 return getP();
572 };

◆ getMomentumVertexErrorMatrix()

TMatrixFSym getMomentumVertexErrorMatrix ( ) const

Returns 7x7 error matrix.

Returns
7x7 error matrix (order: px,py,pz,E,x,y,z)

Definition at line 420 of file Particle.cc.

421{
422 TMatrixFSym m(c_DimMatrix);
423
424 int element = 0;
425 for (int irow = 0; irow < c_DimMatrix; irow++) {
426 for (int icol = irow; icol < c_DimMatrix; icol++) {
427 m(irow, icol) = m(icol, irow) = m_errMatrix[element];
428 element++;
429 }
430 }
431 return m;
432}
double m_errMatrix[c_SizeMatrix]
error matrix (1D representation)
Definition: Particle.h:1067

◆ getMostLikelyTrackFitResult()

std::pair< Const::ChargedStable, const TrackFitResult * > getMostLikelyTrackFitResult ( ) const

For a (track-based) particle, returns the charged stable mass hypothesis associated to the most probable TrackFitResult, and the TrackFitResult itself.

Definition at line 1398 of file Particle.cc.

1399{
1400
1401 const auto track = this->getTrack();
1402
1403 if (!track) {
1404 return std::make_pair(Const::ChargedStable(std::abs(this->getPDGCode())), nullptr);
1405 }
1406
1407 // Find the track fit with the highest pValue
1408 auto trackFitResults = track->getTrackFitResults();
1409 auto it_maxPValue = std::max_element(std::begin(trackFitResults), std::end(trackFitResults),
1410 [](auto tfrAndM1, auto tfrAndM2) {return tfrAndM1.second->getPValue() < tfrAndM2.second->getPValue();});
1411
1412 return trackFitResults[std::distance(std::begin(trackFitResults), it_maxPValue)];
1413
1414}
Provides a type-safe way to pass members of the chargedStableSet set.
Definition: Const.h:589
const Track * getTrack() const
Returns the pointer to the Track object that was used to create this Particle (ParticleType == c_Trac...
Definition: Particle.cc:845

◆ getName()

std::string getName ( ) const
overridevirtual

Return name of this particle.

Reimplemented from RelationsInterface< BASE >.

Definition at line 1165 of file Particle.cc.

1166{
1167 return TDatabasePDG::Instance()->GetParticle(m_pdgCode)->GetName();
1168}

◆ getNDaughters()

unsigned getNDaughters ( void  ) const
inline

Returns number of daughter particles.

Returns
number of daughter particles

Definition at line 727 of file Particle.h.

728 {
729 return m_daughterIndices.size();
730 }

◆ getP()

double getP ( ) const
inline

Returns momentum magnitude (same as getMomentumMagnitude but with shorter name)

Returns
momentum magnitude

Definition at line 578 of file Particle.h.

579 {
581 };

◆ getParticleFromGeneralizedIndexString()

const Particle * getParticleFromGeneralizedIndexString ( const std::string &  generalizedIndex) const

Explores the decay tree of the particle and returns the (grand^n)daughter identified by a generalized index.

The generalized index consists of a colon-separated list of daughter indexes, starting from the root particle: 0:1:3 identifies the fourth daughter (3) of the second daughter (1) of the first daughter (0) of the mother particle.

Parameters
generalizedIndexthe generalized index of the particle to be returned
Returns
a particle in the decay tree of the root particle.

Definition at line 956 of file Particle.cc.

957{
958 // Split the generalizedIndex string in a vector of strings.
959 std::vector<std::string> generalizedIndexes;
960 boost::split(generalizedIndexes, generalizedIndex, boost::is_any_of(":"));
961
962 if (generalizedIndexes.empty()) {
963 B2WARNING("Generalized index of daughter particle is empty. Skipping.");
964 return nullptr;
965 }
966
967 // To explore a decay tree of unknown depth, we need a place to store
968 // both the root particle and the daughter particle at each iteration
969 const Particle* dauPart =
970 nullptr; // This will be eventually returned
971 const Particle* currentPart = this; // This is the root particle of the next iteration
972
973 // Loop over the generalizedIndexes until you get to the particle you want
974 for (auto& indexString : generalizedIndexes) {
975 // indexString is a string. First try to convert it into an int
976 int dauIndex = 0;
977 try {
978 dauIndex = Belle2::convertString<int>(indexString);
979 } catch (std::invalid_argument&) {
980 B2WARNING("Found the string " << indexString << "instead of a daughter index.");
981 return nullptr;
982 }
983
984 // Check that the daughter index is smaller than the number of daughters of the current root particle
985 if (dauIndex >= int(currentPart->getNDaughters()) or dauIndex < 0) {
986 B2WARNING("Daughter index out of range" << LogVar("Daughter index", dauIndex));
987 B2WARNING("Trying to access non-existing particle.");
988 return nullptr;
989 } else {
990 dauPart = currentPart->getDaughter(dauIndex); // Pick the particle indicated by the generalizedIndex
991 currentPart = dauPart;
992 }
993 }
994 return dauPart;
995}

◆ getParticleSource()

EParticleSourceObject getParticleSource ( ) const
inline

Returns particle source as defined with enum EParticleSourceObject.

Returns
particle type

Definition at line 478 of file Particle.h.

479 {
480 return m_particleSource;
481 }

◆ getPDGCode()

int getPDGCode ( void  ) const
inline

Returns PDG code.

Returns
PDG code

Definition at line 454 of file Particle.h.

455 {
456 return m_pdgCode;
457 }

◆ getPDGCodeUsedForFit()

int getPDGCodeUsedForFit ( ) const
inline

Return the always positive PDG code which was used for the track fit (if there was a track fit) of this particle.

This can be different than the Particle's PDG id as not all mass hypothesis are fitted during the reconstruction.

Definition at line 965 of file Particle.h.

966 {
967 return std::abs(m_pdgCodeUsedForFit);
968 }

◆ getPDGLifetime()

double getPDGLifetime ( ) const

Returns particle nominal lifetime.

Returns
nominal lifetime [sec]

Definition at line 613 of file Particle.cc.

614{
615 if (TDatabasePDG::Instance()->GetParticle(m_pdgCode) == nullptr) {
616 B2ERROR("PDG=" << m_pdgCode << " ***code unknown to TDatabasePDG");
617 return 0.0;
618 }
619 return TDatabasePDG::Instance()->GetParticle(m_pdgCode)->Lifetime();
620}

◆ getPDGMass()

double getPDGMass ( void  ) const

Returns uncertainty on the invariant mass (requires valid momentum error matrix)

Returns
mass error Returns particle nominal mass
nominal mass

Definition at line 604 of file Particle.cc.

605{
606 if (TDatabasePDG::Instance()->GetParticle(m_pdgCode) == nullptr) {
607 B2ERROR("PDG=" << m_pdgCode << " ***code unknown to TDatabasePDG");
608 return 0.0;
609 }
610 return TDatabasePDG::Instance()->GetParticle(m_pdgCode)->Mass();
611}

◆ getPIDLikelihood()

const PIDLikelihood * getPIDLikelihood ( ) const

Returns the pointer to the PIDLikelihood object that is related to the Track, which was used to create this Particle (ParticleType == c_Track).

NULL pointer is returned, if the Particle was not made from Track or if the Track has no relation to the PIDLikelihood

Returns
const pointer to the Track

Definition at line 871 of file Particle.cc.

872{
873 if (m_particleSource == c_Track) {
874 StoreArray<Track> tracks{};
875 return tracks[m_mdstIndex]->getRelated<PIDLikelihood>();
876 } else
877 return nullptr;
878}
Class to collect log likelihoods from TOP, ARICH, dEdx, ECL and KLM aimed for output to mdst includes...
Definition: PIDLikelihood.h:29

◆ getProperty()

int getProperty ( ) const
inline

Returns particle property as a bit pattern The values are defined in the PropertyFlags enum and described in detail there.

Returns
Combination of Properties describing the particle property

Definition at line 498 of file Particle.h.

499 {
500 return m_properties;
501 }

◆ getPValue()

double getPValue ( ) const
inline

Returns chi^2 probability of fit if done or -1.

Returns
p-value of fit (nan means no fit done)

Definition at line 667 of file Particle.h.

668 {
669 return m_pValue;
670 }

◆ getPx()

double getPx ( ) const
inline

Returns x component of momentum.

Returns
x component of momentum

Definition at line 587 of file Particle.h.

588 {
590 }

◆ getPy()

double getPy ( ) const
inline

Returns y component of momentum.

Returns
y component of momentum

Definition at line 596 of file Particle.h.

597 {
599 }

◆ getPz()

double getPz ( ) const
inline

Returns z component of momentum.

Returns
z component of momentum

Definition at line 605 of file Particle.h.

606 {
608 }

◆ getRelated()

T * getRelated ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get the object to or from which this object has a relation.

Template Parameters
TThe class of objects to or from which the relation points.
Parameters
nameThe name of the store array to or from which the relation points. If empty the default store array name for class T will be used. If the special name "ALL" is given all store arrays containing objects of type T are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
The first related object or a null pointer.

Definition at line 278 of file RelationsObject.h.

279 {
281 T::Class(), name, namedRelation).object);
282 }
@ c_BothSides
Combination of c_FromSide and c_ToSide.
Definition: DataStore.h:79
Belle2::RelationEntry getRelationWith(ESearchSide searchSide, const TObject *object, StoreEntry *&entry, int &index, const TClass *withClass, const std::string &withName, const std::string &namedRelation)
Get the first relation between an object and another object in a store array.
Definition: DataStore.cc:597
TObject * object
Pointer to the object.
Definition: RelationEntry.h:32

◆ getRelatedFrom()

FROM * getRelatedFrom ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get the object from which this object has a relation.

Template Parameters
FROMThe class of objects from which the relation points.
Parameters
nameThe name of the store array from which the relation points. If empty the default store array name for class FROM will be used. If the special name "ALL" is given all store arrays containing objects of type FROM are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
The first related object or a null pointer.

Definition at line 263 of file RelationsObject.h.

264 {
266 m_cacheArrayIndex, FROM::Class(), name, namedRelation).object);
267 }
@ c_FromSide
Return relations/objects pointed from (to a given object).
Definition: DataStore.h:77

◆ getRelatedFromWithWeight()

std::pair< FROM *, float > getRelatedFromWithWeight ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get first related object & weight of relation pointing from an array.

Template Parameters
FROMThe class of objects from which the relation points.
Parameters
nameThe name of the store array from which the relation points. If empty the default store array name for class FROM will be used. If the special name "ALL" is given all store arrays containing objects of type FROM are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
Pair of first related object and the relation weight, or (NULL, 1.0) if none found.

Definition at line 314 of file RelationsObject.h.

316 {
318 FROM::Class(), name, namedRelation);
319 return std::make_pair(static_cast<FROM*>(entry.object), entry.weight);
320 }

◆ getRelatedTo()

TO * getRelatedTo ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get the object to which this object has a relation.

Template Parameters
TOThe class of objects to which the relation points.
Parameters
nameThe name of the store array to which the relation points. If empty the default store array name for class TO will be used. If the special name "ALL" is given all store arrays containing objects of type TO are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
The first related object or a null pointer.

Definition at line 248 of file RelationsObject.h.

249 {
251 TO::Class(), name, namedRelation).object);
252 }
@ c_ToSide
Return relations/objects pointed to (from a given object).
Definition: DataStore.h:78

◆ getRelatedToWithWeight()

std::pair< TO *, float > getRelatedToWithWeight ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get first related object & weight of relation pointing to an array.

Template Parameters
TOThe class of objects to which the relation points.
Parameters
nameThe name of the store array to which the relation points. If empty the default store array name for class TO will be used. If the special name "ALL" is given all store arrays containing objects of type TO are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
Pair of first related object and the relation weight, or (NULL, 1.0) if none found.

Definition at line 297 of file RelationsObject.h.

299 {
301 TO::Class(), name, namedRelation);
302 return std::make_pair(static_cast<TO*>(entry.object), entry.weight);
303 }

◆ getRelatedWithWeight()

std::pair< T *, float > getRelatedWithWeight ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get first related object & weight of relation pointing from/to an array.

Template Parameters
TThe class of objects to or from which the relation points.
Parameters
nameThe name of the store array to or from which the relation points. If empty the default store array name for class T will be used. If the special name "ALL" is given all store arrays containing objects of type T are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
Pair of first related object and the relation weight, or (NULL, 1.0) if none found.

Definition at line 331 of file RelationsObject.h.

333 {
335 T::Class(), name, namedRelation);
336 return std::make_pair(static_cast<T*>(entry.object), entry.weight);
337 }

◆ getRelationsFrom()

RelationVector< FROM > getRelationsFrom ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get the relations that point from another store array to this object.

Template Parameters
FROMThe class of objects from which the relations point.
Parameters
nameThe name of the store array from which the relations point. If empty the default store array name for class FROM will be used. If the special name "ALL" is given all store arrays containing objects of type FROM are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
A vector of relations.

Definition at line 212 of file RelationsObject.h.

214 {
216 m_cacheArrayIndex, FROM::Class(), name, namedRelation));
217 }
RelationVector< T > getRelationsWith(const std::string &name="", const std::string &namedRelation="") const
Get the relations between this object and another store array.

◆ getRelationsTo()

RelationVector< TO > getRelationsTo ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get the relations that point from this object to another store array.

Template Parameters
TOThe class of objects to which the relations point.
Parameters
nameThe name of the store array to which the relations point. If empty the default store array name for class TO will be used. If the special name "ALL" is given all store arrays containing objects of type TO are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
A vector of relations.

Definition at line 197 of file RelationsObject.h.

198 {
200 m_cacheArrayIndex, TO::Class(), name, namedRelation));
201 }

◆ getRelationsWith()

RelationVector< T > getRelationsWith ( const std::string &  name = "",
const std::string &  namedRelation = "" 
) const
inlineinherited

Get the relations between this object and another store array.

Relations in both directions are returned.

Template Parameters
TThe class of objects to or from which the relations point.
Parameters
nameThe name of the store array to or from which the relations point. If empty the default store array name for class T will be used. If the special name "ALL" is given all store arrays containing objects of type T are considered.
namedRelationAdditional name for the relation, or "" for the default naming
Returns
A vector of relations.

Definition at line 230 of file RelationsObject.h.

231 {
233 m_cacheArrayIndex, T::Class(), name, namedRelation));
234 }

◆ getTrack()

const Track * getTrack ( ) const

Returns the pointer to the Track object that was used to create this Particle (ParticleType == c_Track).

NULL pointer is returned, if the Particle was not made from Track.

Returns
const pointer to the Track

Definition at line 845 of file Particle.cc.

846{
847 if (m_particleSource == c_Track) {
848 StoreArray<Track> tracks{};
849 return tracks[m_mdstIndex];
850 } else
851 return nullptr;
852}

◆ getTrackFitResult()

const TrackFitResult * getTrackFitResult ( ) const

Returns the pointer to the TrackFitResult that was used to create this Particle (ParticleType == c_Track).

NULL pointer is returned, if the Particle was not made from Track.

Returns
const pointer to the TrackFitResult

Definition at line 854 of file Particle.cc.

855{
856 // if the particle is related to a TrackFitResult then return this
857 auto* selfrelated = this->getRelatedTo<TrackFitResult>();
858 if (selfrelated && !isnan(selfrelated->getPValue()))
859 return selfrelated;
860
861 // if not get the TFR with closest mass to this particle
862 auto* selftrack = this->getTrack();
863 if (selftrack)
864 return selftrack->getTrackFitResultWithClosestMass(
865 Belle2::Const::ChargedStable(std::abs(this->getPDGCode())));
866
867 // otherwise we're probably not a track based particle
868 return nullptr;
869}

◆ getV0()

const V0 * getV0 ( ) const

Returns the pointer to the V0 object that was used to create this Particle (if ParticleType == c_V0).

NULL pointer is returned if the Particle was not made from a V0.

Returns
const pointer to the V0

Definition at line 880 of file Particle.cc.

881{
882 if (m_particleSource == c_V0) {
883 StoreArray<V0> v0s{};
884 return v0s[m_mdstIndex];
885 } else {
886 return nullptr;
887 }
888}

◆ getVertex()

ROOT::Math::XYZVector getVertex ( ) const
inline

Returns vertex position (POCA for charged, IP for neutral FS particles)

Returns
vertex position

Definition at line 631 of file Particle.h.

632 {
633 return ROOT::Math::XYZVector(m_x, m_y, m_z);
634 };

◆ getVertexErrorMatrix()

TMatrixFSym getVertexErrorMatrix ( ) const

Returns the 3x3 position error sub-matrix.

Returns
The 3x3 position error matrix (order: x,y,z)

Definition at line 447 of file Particle.cc.

448{
449 TMatrixFSym pos;
450 const TMatrixFSym& full = getMomentumVertexErrorMatrix();
451
452 // get 3x3 (position) submatrix from the full error matrix
453 // vertex related elements are in [4,5,6]x[4,5,6] block
454 full.GetSub(4, 6, pos, "S");
455
456 return pos;
457}

◆ getX()

double getX ( ) const
inline

Returns x component of vertex position.

Returns
x component of vertex position

Definition at line 640 of file Particle.h.

641 {
642 return m_x;
643 }

◆ getY()

double getY ( ) const
inline

Returns y component of vertex position.

Returns
y component of vertex position

Definition at line 649 of file Particle.h.

650 {
651 return m_y;
652 }

◆ getZ()

double getZ ( ) const
inline

Returns z component of vertex position.

Returns
z component of vertex position

Definition at line 658 of file Particle.h.

659 {
660 return m_z;
661 }

◆ hasExtraInfo()

bool hasExtraInfo ( const std::string &  name) const

Return whether the extra info with the given name is set.

Definition at line 1266 of file Particle.cc.

1267{
1268 if (m_extraInfo.empty())
1269 return false;
1270
1271 //get index for name
1272 const auto mapID = (unsigned int)m_extraInfo[0];
1274 if (!extraInfoMap) {
1275 B2FATAL("ParticleExtraInfoMap not available, but needed for storing extra info in Particle!");
1276 }
1277 unsigned int index = extraInfoMap->getIndex(mapID, name);
1278 if (index == 0 or index >= m_extraInfo.size()) //actually indices start at 1
1279 return false;
1280
1281 return true;
1282}

◆ isCopyOf()

bool isCopyOf ( const Particle oParticle,
bool  doDetailedComparison = false 
) const

Returns true if this Particle and oParticle are copies of each other.

Copies are defined as: o) The decay chain of both Particles is exactly the same o) Both Particles are constructed from exactly the same final state particles

  • Examples: M1 -> (C1 -> F1 F2) (C2 -> F3 F4) M2 -> (C1 -> F1 F2) (C2 -> F3 F5) M3 -> (C1 -> F1 F2) F3 F4 M4 -> (C1 -> F1 F2) (C2 -> F3 F4) and M4 is kinematically fitted (its momentum updated)

M1 and M2 are not copies since condition 2 is not fulfilled. M1 and M3 are not copies since condition 1 is not fulfilled. M1 and M4 are copies since both conditions are fulfilled.

Parameters
oParticlepointer to other particle
doDetailedComparisonif true, this means that particles of different PDG codes, but created from the same track or cluster will be indicated as copies. Returns B2FATAL in case of comparison of c_MCParticle type to a non c_MCParticle.
Returns
true if particles are copies of each-other, otherwise false

Definition at line 752 of file Particle.cc.

753{
754 // the name of the game is to as quickly as possible determine
755 // that the Particles are not copies
756 if (this->getPDGCode() != oParticle->getPDGCode() and !doDetailedComparison)
757 return false;
758
759 unsigned nDaughters = this->getNDaughters();
760 if (nDaughters != oParticle->getNDaughters())
761 return false;
762
763 if (nDaughters) {
764 // has daughters: check if the decay chain is the same and it ends with
765 // the same FSPs
766 std::vector<int> thisDecayChain(nDaughters * 2);
767 std::vector<int> othrDecayChain(nDaughters * 2);
768
769 for (unsigned i = 0; i < nDaughters; i++) {
770 this->getDaughter(i)->fillDecayChain(thisDecayChain);
771 oParticle->getDaughter(i)->fillDecayChain(othrDecayChain);
772 }
773
774 // Even if the daughters have been provided in a different order in the
775 // decay string the particles should be identified as copies / duplicates.
776 // Therefore, the decay chain vectors, which contain both the pdg codes and
777 // the mdst sources, are sorted here before comparing them.
778 sort(thisDecayChain.begin(), thisDecayChain.end());
779 sort(othrDecayChain.begin(), othrDecayChain.end());
780
781 for (unsigned i = 0; i < thisDecayChain.size(); i++)
782 if (thisDecayChain[i] != othrDecayChain[i])
783 return false;
784
785 } else if (this->getMdstSource() != oParticle->getMdstSource() and !doDetailedComparison) {
786 // has no daughters: it's a FSP, compare MDST source and index
787 return false;
788 }
789 // Stop here if we do not want a detailed comparison
790 if (!doDetailedComparison)
791 return true;
792 //If we compare here a reconstructed Particle to a generated MCParticle
793 //it means that something went horribly wrong and we must stop
794 if ((this->getParticleSource() == EParticleSourceObject::c_MCParticle
795 and oParticle->getParticleSource() != EParticleSourceObject::c_MCParticle)
796 or (this->getParticleSource() != EParticleSourceObject::c_MCParticle
797 and oParticle->getParticleSource() == EParticleSourceObject::c_MCParticle)) {
798 B2WARNING("Something went wrong: MCParticle is being compared to a non MC Particle. Please check your script!\n"
799 " If the MCParticle <-> Particle comparison happens in the RestOfEventBuilder,\n"
800 " the Rest Of Event may contain signal side particles.");
801 return false;
802 }
803 if (this->getParticleSource() == EParticleSourceObject::c_MCParticle
804 and oParticle->getParticleSource() == EParticleSourceObject::c_MCParticle) {
805 return this->getMCParticle() == oParticle->getMCParticle();
806 }
807 if (this->getParticleSource() != oParticle->getParticleSource()) {
808 return false;
809 }
810 if (this->getMdstSource() == oParticle->getMdstSource()) {
811 return true;
812 }
813 if (this->getTrack() && oParticle->getTrack() &&
814 this->getTrack()->getArrayIndex() != oParticle->getTrack()->getArrayIndex()) {
815 return false;
816 }
817 if (this->getKLMCluster() && oParticle->getKLMCluster()
818 && this->getKLMCluster()->getArrayIndex() != oParticle->getKLMCluster()->getArrayIndex()) {
819 return false;
820 }
821
822 // It can be a bit more complicated for ECLClusters as we might also have to ensure they are connected-region unique
823 if (this->getECLCluster() && oParticle->getECLCluster()
824 && this->getECLCluster()->getArrayIndex() != oParticle->getECLCluster()->getArrayIndex()) {
825
826 // if either is a track then they must be different
827 if (this->getECLCluster()->isTrack() or oParticle->getECLCluster()->isTrack())
828 return false;
829
830 // we cannot combine two particles of different hypotheses from the same
831 // connected region (as their energies overlap)
833 return false;
834
835 // in the rare case that both are neutral and the hypotheses are different,
836 // we must also check that they are from different connected regions
837 // otherwise they come from the "same" underlying ECLShower
838 if (this->getECLCluster()->getConnectedRegionId() != oParticle->getECLCluster()->getConnectedRegionId())
839 return false;
840 }
841 return true;
842
843}
bool isTrack() const
Return true if the cluster matches with track.
Definition: ECLCluster.h:235
int getConnectedRegionId() const
Return connected region id.
Definition: ECLCluster.h:250
const KLMCluster * getKLMCluster() const
Returns the pointer to the KLMCluster object that was used to create this Particle (ParticleType == c...
Definition: Particle.cc:926
const MCParticle * getMCParticle() const
Returns the pointer to the MCParticle object that was used to create this Particle (ParticleType == c...
Definition: Particle.cc:942
EParticleSourceObject getParticleSource() const
Returns particle source as defined with enum EParticleSourceObject.
Definition: Particle.h:478

◆ isMostLikely()

bool isMostLikely ( ) const

Returns true if the (track-based) particle is created with its most likely mass hypothesis based on PID likelihood.

Definition at line 1391 of file Particle.cc.

1392{
1393 const PIDLikelihood* likelihood = Particle::getPIDLikelihood();
1394 if (likelihood) return likelihood->getMostLikely().getPDGCode() == std::abs(m_pdgCode);
1395 else return false;
1396}
Const::ChargedStable getMostLikely(const double *fractions=nullptr, Const::PIDDetectorSet set=Const::PIDDetectorSet::set()) const
Return most likely particle among chargedStableSet; if prior fractions not given equal prior probabil...
const PIDLikelihood * getPIDLikelihood() const
Returns the pointer to the PIDLikelihood object that is related to the Track, which was used to creat...
Definition: Particle.cc:871

◆ isMostLikelyTrackFitResult()

bool isMostLikelyTrackFitResult ( ) const

Returns true if the (track-based) particle is created with its most likely mass hypothesis based on TrackFitResult.

Definition at line 1416 of file Particle.cc.

1417{
1418 const auto trackFit = this->getTrackFitResult();
1419 if (!trackFit) return false;
1420
1421 return (trackFit == this->getMostLikelyTrackFitResult().second);
1422}
std::pair< Const::ChargedStable, const TrackFitResult * > getMostLikelyTrackFitResult() const
For a (track-based) particle, returns the charged stable mass hypothesis associated to the most proba...
Definition: Particle.cc:1398
const TrackFitResult * getTrackFitResult() const
Returns the pointer to the TrackFitResult that was used to create this Particle (ParticleType == c_Tr...
Definition: Particle.cc:854

◆ overlapsWith()

bool overlapsWith ( const Particle oParticle) const

Returns true if final state ancestors of oParticle overlap.

Parameters
oParticlepointer to particle
Returns
true if overlap, otherwise false

Definition at line 737 of file Particle.cc.

738{
739 // obtain vectors of daughter final state particles
740 std::vector<const Particle*> thisFSPs = this->getFinalStateDaughters();
741 std::vector<const Particle*> otherFSPs = oParticle->getFinalStateDaughters();
742
743 // check if they share any of the FSPs
744 for (auto& thisFSP : thisFSPs)
745 for (auto& otherFSP : otherFSPs)
746 if (thisFSP->getMdstSource() == otherFSP->getMdstSource())
747 return true;
748
749 return false;
750}

◆ print()

void print ( ) const

Prints the contents of a Particle object to standard output.

Definition at line 1170 of file Particle.cc.

1171{
1172 B2INFO(getInfo());
1173}
std::string getInfo() const
Return a short summary of this object's contents in raw text format.

◆ removeDaughter()

void removeDaughter ( const Particle daughter,
const bool  updateType = true 
)

Removes index of daughter from daughters index array.

Parameters
daughterpointer to the daughter particle
updateTypebool whether particle type should be updated if last daughter was removed

Definition at line 688 of file Particle.cc.

689{
690 if (getNDaughters() == 0)
691 return;
692
693 for (unsigned i = 0; i < getNDaughters(); i++) {
694 if (m_daughterIndices[i] == daughter->getArrayIndex()) {
695 m_daughterIndices.erase(m_daughterIndices.begin() + i);
697 i--;
698 }
699 }
700
701 if (getNDaughters() == 0 and updateType)
702 m_particleSource = c_Undefined;
703}

◆ removeExtraInfo()

void removeExtraInfo ( )

Remove all stored extra info fields.

Definition at line 1284 of file Particle.cc.

1285{
1286 m_extraInfo.clear();
1287}

◆ replaceDaughter()

bool replaceDaughter ( const Particle oldDaughter,
Particle newDaughter 
)

Replace index of given daughter with new daughter, return true if a replacement is made.

Parameters
oldDaughterpointer to the daughter that will be removed
newDaughterpointer to the particle that will be added as a daughter

Definition at line 705 of file Particle.cc.

706{
707 int index = oldDaughter->getArrayIndex();
708
709 for (unsigned i = 0; i < getNDaughters(); i++) {
710 if (m_daughterIndices[i] == index) {
711 auto ite_index = m_daughterIndices.erase(m_daughterIndices.begin() + i);
712 m_daughterIndices.insert(ite_index, newDaughter->getArrayIndex());
713 auto ite_property = m_daughterProperties.erase(m_daughterProperties.begin() + i);
714 m_daughterProperties.insert(ite_property, Particle::PropertyFlags::c_Ordinary);
715
716 newDaughter->writeExtraInfo("original_index", index);
717 return true;
718 }
719 }
720 return false;
721}
void writeExtraInfo(const std::string &name, const double value)
Sets the user defined extraInfo.
Definition: Particle.cc:1308

◆ replaceDaughterRecursively()

bool replaceDaughterRecursively ( const Particle oldDaughter,
Particle newDaughter 
)

Apply replaceDaughter to all Particles in the decay tree by looping recursively through it, return true if a replacement is made.

Parameters
oldDaughterpointer to the daughter that will be removed
newDaughterpointer to the particle that will be added as a daughter

Definition at line 723 of file Particle.cc.

724{
725 bool isReplaced = this->replaceDaughter(oldDaughter, newDaughter);
726 if (isReplaced)
727 return true;
728 for (auto& daughter : this->getDaughters()) {
729 isReplaced = daughter->replaceDaughterRecursively(oldDaughter, newDaughter);
730 if (isReplaced)
731 return true;
732 }
733 return false;
734}
bool replaceDaughter(const Particle *oldDaughter, Particle *newDaughter)
Replace index of given daughter with new daughter, return true if a replacement is made.
Definition: Particle.cc:705

◆ resetErrorMatrix()

void resetErrorMatrix ( )
private

Resets 7x7 error matrix All elements are set to 0.0.

Definition at line 1079 of file Particle.cc.

1080{
1081 for (double& i : m_errMatrix)
1082 i = 0.0;
1083}

◆ resetJacobiMatrix()

void resetJacobiMatrix ( )
private

Resets 4x6 error matrix All elements are set to 0.0.

Definition at line 1085 of file Particle.cc.

1086{
1087 for (double& i : m_jacobiMatrix)
1088 i = 0.0;
1089}
double m_jacobiMatrix[c_SizeMatrix]
error matrix (1D representation)
Definition: Particle.h:1068

◆ set4Vector()

void set4Vector ( const ROOT::Math::PxPyPzEVector &  p4)
inline

Sets Lorentz vector.

Parameters
p4Lorentz vector

Definition at line 271 of file Particle.h.

272 {
273 m_px = p4.Px();
274 m_py = p4.Py();
275 m_pz = p4.Pz();
276 m_mass = p4.M();
277 }

◆ set4VectorDividingByMomentumScaling()

void set4VectorDividingByMomentumScaling ( const ROOT::Math::PxPyPzEVector &  p4)
inline

Sets Lorentz vector dividing by the momentum scaling factor.

Parameters
p4Lorentz vector

Definition at line 283 of file Particle.h.

284 {
285 m_px = p4.Px() / m_momentumScale;
286 m_py = p4.Py() / m_momentumScale;
287 m_pz = p4.Pz() / m_momentumScale;
288 m_mass = p4.M();
289 }

◆ setEnergyLossCorrection()

void setEnergyLossCorrection ( double  energyLossCorrection)
inline

Sets Energy loss correction.

Parameters
energyLossCorrectionCorrection factor

Definition at line 306 of file Particle.h.

307 {
308 m_energyLossCorrection = energyLossCorrection;
309 }

◆ setExtraInfo()

void setExtraInfo ( const std::string &  name,
double  value 
)

Sets the user-defined data of given name to the given value.

throws std::runtime_error if variable isn't set.

Definition at line 1317 of file Particle.cc.

1318{
1319 if (m_extraInfo.empty())
1320 throw std::runtime_error(std::string("setExtraInfo: Value '") + name + "' not found in Particle!");
1321
1322 //get index for name
1323 const auto mapID = (unsigned int)m_extraInfo[0];
1325 if (!extraInfoMap) {
1326 B2FATAL("ParticleExtraInfoMap not available, but needed for storing extra info in Particle!");
1327 }
1328 unsigned int index = extraInfoMap->getIndex(mapID, name);
1329 if (index == 0 or index >= m_extraInfo.size()) //actually indices start at 1
1330 throw std::runtime_error(std::string("setExtraInfo: Value '") + name + "' not found in Particle!");
1331
1332 m_extraInfo[index] = value;
1333
1334}

◆ setFlavorType()

void setFlavorType ( )
private

sets m_flavorType using m_pdgCode

Definition at line 1150 of file Particle.cc.

1151{
1153 if (m_pdgCode < 0) return;
1154 if (m_pdgCode == Const::photon.getPDGCode()) {m_flavorType = c_Unflavored; return;} // gamma
1155 if (m_pdgCode == Const::Kshort.getPDGCode()) {m_flavorType = c_Unflavored; return;} // K_s
1156 if (m_pdgCode == Const::Klong.getPDGCode()) {m_flavorType = c_Unflavored; return;} // K_L
1157 if (m_pdgCode == 43) {m_flavorType = c_Unflavored; return;} // Xu0
1158 int nnn = m_pdgCode / 10;
1159 int q3 = nnn % 10; nnn /= 10;
1160 int q2 = nnn % 10; nnn /= 10;
1161 int q1 = nnn % 10;
1162 if (q1 == 0 && q2 == q3) m_flavorType = c_Unflavored; // unflavored meson
1163}
static const ParticleType Kshort
K^0_S particle.
Definition: Const.h:677

◆ setJacobiMatrix()

void setJacobiMatrix ( const TMatrixF &  jacobiMatrix)

Sets 4x6 jacobi matrix.

Parameters
jacobiMatrix4x6 momentum and vertex error matrix (order: px,py,pz,E,x,y,z)

Definition at line 407 of file Particle.cc.

408{
409 // check if provided Jacobi Matrix is of dimension 4x6
410 // if not, reset the error matrix and print warning
411 if (m.GetNrows() != 4 || m.GetNcols() != 6) {
413 B2WARNING("Jacobi Matrix is not 4x6 ");
414 return;
415 }
417}
void resetJacobiMatrix()
Resets 4x6 error matrix All elements are set to 0.0.
Definition: Particle.cc:1085
void storeJacobiMatrix(const TMatrixF &jacobiMatrix)
Stores 4x6 Jacobi matrix into private member m_jacobiMatrix.
Definition: Particle.cc:1102

◆ setMdstArrayIndex()

void setMdstArrayIndex ( const int  arrayIndex)
private

set mdst array index

Definition at line 350 of file Particle.cc.

351{
352 m_mdstIndex = arrayIndex;
353
354 // set the identifier
355 if (m_particleSource == c_ECLCluster) {
356 const ECLCluster* cluster = this->getECLCluster();
357 if (cluster) {
358 m_identifier = cluster->getUniqueClusterId();
359 } else {
360 B2ERROR("Particle is of type = ECLCluster has identifier not set and no relation to ECLCluster.\n"
361 "This has happen because old microDST is analysed with newer version of software.");
362 }
363 } else {
365 }
366}

◆ setMomentumPositionErrorMatrix()

void setMomentumPositionErrorMatrix ( const TrackFitResult trackFit)
private

Sets the momentum, position and error matrix for this particle (created from charged Track)

Definition at line 1001 of file Particle.cc.

1002{
1003 // set momentum
1004 m_px = trackFit->getMomentum().X();
1005 m_py = trackFit->getMomentum().Y();
1006 m_pz = trackFit->getMomentum().Z();
1007
1008 // set position at which the momentum is given (= POCA)
1009 setVertex(trackFit->getPosition());
1010
1011 // set Chi^2 probability
1012 m_pValue = trackFit->getPValue();
1013
1014 // set error matrix
1015 const auto cov6 = trackFit->getCovariance6();
1016 constexpr unsigned order[] = {c_X, c_Y, c_Z, c_Px, c_Py, c_Pz};
1017
1018 TMatrixFSym errMatrix(c_DimMatrix);
1019 for (int i = 0; i < 6; i++) {
1020 for (int j = i; j < 6; j++) {
1021 // although it seems to make no sense to fill all elements of the
1022 // symmetric matrix, it has to be (do not touch this code)
1023 errMatrix(order[j], order[i]) = errMatrix(order[i], order[j]) = cov6(i, j);
1024 }
1025 }
1026
1027 /*
1028 E = sqrt(px^2 + py^2 + pz^2 + m^2) thus:
1029 cov(x,E) = cov(px,x) *dE/dpx + cov(py,x) *dE/dpy + cov(pz,x) *dE/dpz
1030 cov(y,E) = cov(px,y) *dE/dpx + cov(py,y) *dE/dpy + cov(pz,y) *dE/dpz
1031 cov(z,E) = cov(px,z) *dE/dpx + cov(py,z) *dE/dpy + cov(pz,z) *dE/dpz
1032 cov(px,E) = cov(px,px)*dE/dpx + cov(px,py)*dE/dpy + cov(px,pz)*dE/dpz
1033 cov(py,E) = cov(py,px)*dE/dpx + cov(py,py)*dE/dpy + cov(py,pz)*dE/dpz
1034 cov(pz,E) = cov(pz,px)*dE/dpx + cov(pz,py)*dE/dpy + cov(pz,pz)*dE/dpz
1035 cov(E,E) = cov(px,px)*(dE/dpx)^2 + cov(py,py)*(dE/dpy)^2 + cov(pz,pz)*(dE/dpz)^2
1036 + 2*cov(px,py)*dE/dpx*dE/dpy
1037 + 2*cov(py,pz)*dE/dpy*dE/dpz
1038 + 2*cov(pz,px)*dE/dpz*dE/dpx
1039 dE/dpx = px/E etc.
1040 */
1041
1042 const double E = getEnergy();
1043 const double dEdp[] = {m_px / E, m_py / E, m_pz / E};
1044 constexpr unsigned compMom[] = {c_Px, c_Py, c_Pz};
1045 constexpr unsigned compPos[] = {c_X, c_Y, c_Z};
1046
1047 // covariances (p,E)
1048 for (unsigned int i : compMom) {
1049 double Cov = 0;
1050 for (int k = 0; k < 3; k++) {
1051 Cov += errMatrix(i, compMom[k]) * dEdp[k];
1052 }
1053 errMatrix(i, c_E) = Cov;
1054 }
1055
1056 // covariances (x,E)
1057 for (unsigned int comp : compPos) {
1058 double Cov = 0;
1059 for (int k = 0; k < 3; k++) {
1060 Cov += errMatrix(comp, compMom[k]) * dEdp[k];
1061 }
1062 errMatrix(c_E, comp) = Cov;
1063 }
1064
1065 // variance (E,E)
1066 double Cov = 0;
1067 for (int i = 0; i < 3; i++) {
1068 Cov += errMatrix(compMom[i], compMom[i]) * dEdp[i] * dEdp[i];
1069 }
1070 for (int i = 0; i < 3; i++) {
1071 int k = (i + 1) % 3;
1072 Cov += 2 * errMatrix(compMom[i], compMom[k]) * dEdp[i] * dEdp[k];
1073 }
1074 errMatrix(c_E, c_E) = Cov;
1075
1076 storeErrorMatrix(errMatrix);
1077}
R E
internal precision of FFTW codelets
void storeErrorMatrix(const TMatrixFSym &errMatrix)
Stores 7x7 error matrix into private member m_errMatrix.
Definition: Particle.cc:1091
double getPValue() const
Getter for Chi2 Probability of the track fit.
TMatrixDSym getCovariance6() const
Position and Momentum Covariance Matrix.
ROOT::Math::XYZVector getMomentum() const
Getter for vector of momentum at closest approach of track in r/phi projection.
ROOT::Math::XYZVector getPosition() const
Getter for vector of position at closest approach of track in r/phi projection.

◆ setMomentumScalingFactor()

void setMomentumScalingFactor ( double  momentumScalingFactor)
inline

Sets momentum scaling.

Parameters
momentumScalingFactorscaling factor

Definition at line 334 of file Particle.h.

335 {
336 m_momentumScalingFactor = momentumScalingFactor;
338 }
double m_momentumSmearingFactor
momentum smearing factor
Definition: Particle.h:1062
double m_momentumScalingFactor
momentum scaling factor
Definition: Particle.h:1061

◆ setMomentumSmearingFactor()

void setMomentumSmearingFactor ( double  momentumSmearingFactor)
inline

Sets momentum smearing.

Parameters
momentumSmearingFactorscaling factor

Definition at line 344 of file Particle.h.

◆ setMomentumVertexErrorMatrix()

void setMomentumVertexErrorMatrix ( const TMatrixFSym &  errMatrix)

Sets 7x7 error matrix.

Parameters
errMatrix7x7 momentum and vertex error matrix (order: px,py,pz,E,x,y,z)

Definition at line 393 of file Particle.cc.

394{
395 // check if provided Error Matrix is of dimension 7x7
396 // if not, reset the error matrix and print warning
397 if (m.GetNrows() != c_DimMatrix || m.GetNcols() != c_DimMatrix) {
399 B2WARNING("Error Matrix is not 7x7 ");
400 return;
401 }
403
404}

◆ setPDGCode()

void setPDGCode ( const int  pdg)
inline

Sets PDG code.

Parameters
pdgPDG code

Definition at line 262 of file Particle.h.

263 {
264 m_pdgCode = pdg;
265 }

◆ setProperty()

void setProperty ( const int  properties)
inline

sets m_properties

Definition at line 374 of file Particle.h.

375 {
376 m_properties = properties;
377 }

◆ setPValue()

void setPValue ( double  pValue)
inline

Sets chi^2 probability of fit.

Parameters
pValuep-value of fit

Definition at line 366 of file Particle.h.

367 {
368 m_pValue = pValue;
369 }

◆ setVertex()

void setVertex ( const ROOT::Math::XYZVector &  vertex)
inline

Sets position (decay vertex)

Parameters
vertexposition

Definition at line 295 of file Particle.h.

296 {
297 m_x = vertex.X();
298 m_y = vertex.Y();
299 m_z = vertex.Z();
300 };

◆ storeErrorMatrix()

void storeErrorMatrix ( const TMatrixFSym &  errMatrix)
private

Stores 7x7 error matrix into private member m_errMatrix.

Parameters
errMatrix7x7 error matrix

Definition at line 1091 of file Particle.cc.

1092{
1093 int element = 0;
1094 for (int irow = 0; irow < c_DimMatrix; irow++) {
1095 for (int icol = irow; icol < c_DimMatrix; icol++) {
1096 m_errMatrix[element] = m(irow, icol);
1097 element++;
1098 }
1099 }
1100}

◆ storeJacobiMatrix()

void storeJacobiMatrix ( const TMatrixF &  jacobiMatrix)
private

Stores 4x6 Jacobi matrix into private member m_jacobiMatrix.

Parameters
jacobiMatrix4x6 error matrix

Definition at line 1102 of file Particle.cc.

1103{
1104 int element = 0;
1105 for (int irow = 0; irow < 4; irow++) {
1106 for (int icol = 0; icol < 6; icol++) {
1107 m_jacobiMatrix[element] = m(irow, icol);
1108 element++;
1109 }
1110 }
1111}

◆ updateJacobiMatrix()

void updateJacobiMatrix ( )

Propagate the photon energy scaling to jacobian elements that were calculated using energy.

Definition at line 257 of file Particle.cc.

258{
259 ClusterUtils C;
260
261 const ECLCluster* cluster = this->getECLCluster();
262
263 const XYZVector clustervertex = C.GetIPPosition();
264
265 // Get Jacobi matrix.
266 TMatrixD jacobi = C.GetJacobiMatrix4x6FromCluster(cluster, clustervertex, getECLClusterEHypothesisBit());
267 storeJacobiMatrix(jacobi);
268
269 // Propagate the photon energy scaling to jacobian elements that were calculated using energy.
270 TMatrixD scaledJacobi(4, 6);
271
272 int element = 0;
273
274 for (int irow = 0; irow < 4; irow++) {
275 for (int icol = 0; icol < 6; icol++) {
276 if (icol != 0 && irow != 3) {
277 scaledJacobi(irow, icol) = m_jacobiMatrix[element] * m_momentumScale;
278 } else {
279 scaledJacobi(irow, icol) = m_jacobiMatrix[element];
280 }
281 element++;
282 }
283 }
284
285 storeJacobiMatrix(scaledJacobi);
286
287 // Get covariance matrix of IP distribution.
288 const TMatrixDSym clustervertexcovmat = C.GetIPPositionCovarianceMatrix();
289
290 // Set error matrix.
291 TMatrixDSym clustercovmat = C.GetCovarianceMatrix7x7FromCluster(cluster, clustervertexcovmat, scaledJacobi);
292 storeErrorMatrix(clustercovmat);
293}
const TMatrixDSym GetCovarianceMatrix7x7FromCluster(const ECLCluster *cluster, const TMatrixD &jacobiMatrix)
Returns 7x7 covariance matrix (px, py, pz, E, x, y, z)
const TMatrixD GetJacobiMatrix4x6FromCluster(const ECLCluster *cluster, ECLCluster::EHypothesisBit hypo)
Returns 4x6 Jacobi matrix (px, py, pz, E)
Definition: ClusterUtils.cc:54
const TMatrixDSym GetIPPositionCovarianceMatrix()
Returns default IP position covariance matrix from beam parameters.

◆ updateMass()

void updateMass ( const int  pdgCode)

Updates particle mass with the mass of the particle corresponding to the given PDG.

Parameters
pdgCodePDG code of the particle with the desired mass

Definition at line 597 of file Particle.cc.

598{
599 if (TDatabasePDG::Instance()->GetParticle(pdgCode) == nullptr)
600 B2FATAL("PDG=" << pdgCode << " ***code unknown to TDatabasePDG");
601 m_mass = TDatabasePDG::Instance()->GetParticle(pdgCode)->Mass() ;
602}

◆ updateMomentum()

void updateMomentum ( const ROOT::Math::PxPyPzEVector &  p4,
const ROOT::Math::XYZVector &  vertex,
const TMatrixFSym &  errMatrix,
double  pValue 
)
inline

Sets Lorentz vector, position, 7x7 error matrix and p-value.

Parameters
p4Lorentz vector
vertexpoint (position or vertex)
errMatrix7x7 momentum and vertex error matrix (order: px,py,pz,E,x,y,z)
pValuechi^2 probability of the fit

Definition at line 386 of file Particle.h.

390 {
391 set4Vector(p4);
392 setVertex(vertex);
394 m_pValue = pValue;
395 }
void setMomentumVertexErrorMatrix(const TMatrixFSym &errMatrix)
Sets 7x7 error matrix.
Definition: Particle.cc:393

◆ wasExactFitHypothesisUsed()

bool wasExactFitHypothesisUsed ( ) const
inline

Returns true if the type represented by this Particle object was used use as a mass hypothesis during the track of this Particle's parameters.

Definition at line 975 of file Particle.h.

976 {
977 return std::abs(m_pdgCodeUsedForFit) == std::abs(m_pdgCode);
978 }

◆ writeExtraInfo()

void writeExtraInfo ( const std::string &  name,
const double  value 
)

Sets the user defined extraInfo.

Adds it if necessary, overwrites existing ones if they share the same name.

Definition at line 1308 of file Particle.cc.

1309{
1310 if (this->hasExtraInfo(name)) {
1311 this->setExtraInfo(name, value);
1312 } else {
1313 this->addExtraInfo(name, value);
1314 }
1315}
void setExtraInfo(const std::string &name, double value)
Sets the user-defined data of given name to the given value.
Definition: Particle.cc:1317
void addExtraInfo(const std::string &name, double value)
Sets the user-defined data of given name to the given value.
Definition: Particle.cc:1336

Friends And Related Function Documentation

◆ ParticleSubset

friend class ParticleSubset
friend

Definition at line 1170 of file Particle.h.

Member Data Documentation

◆ m_arrayPointer

TClonesArray* m_arrayPointer
mutableprivate

Internal pointer to DataStore array containing the daughters of this particle.

This is a transient member and will not be written to file. The pointer is set by getArrayPointer() when first called.

Definition at line 1100 of file Particle.h.

◆ m_cacheArrayIndex

int m_cacheArrayIndex
mutableprivateinherited

Cache of the index in the TClonesArray to which this object belongs.

Definition at line 432 of file RelationsObject.h.

◆ m_cacheDataStoreEntry

DataStore::StoreEntry* m_cacheDataStoreEntry
mutableprivateinherited

Cache of the data store entry to which this object belongs.

Definition at line 429 of file RelationsObject.h.

◆ m_daughterIndices

std::vector<int> m_daughterIndices
private

daughter particle indices

Definition at line 1070 of file Particle.h.

◆ m_daughterProperties

std::vector<int> m_daughterProperties
private

daughter particle properties

Definition at line 1075 of file Particle.h.

◆ m_energyLossCorrection

double m_energyLossCorrection = 0.0
private

energy loss correction.

defined as 'm_energyLossCorrection = E_measured - E_true'

Definition at line 1063 of file Particle.h.

◆ m_errMatrix

double m_errMatrix[c_SizeMatrix] = {}
private

error matrix (1D representation)

Definition at line 1067 of file Particle.h.

◆ m_extraInfo

std::vector<double> m_extraInfo
private

Stores associated user defined values.

Order is given by string -> index mapping in ParticleExtraInfoMap. entry 0 is reserved specifies which map to use.

Definition at line 1091 of file Particle.h.

◆ m_flavorType

EFlavorType m_flavorType
private

flavor type.

Definition at line 1071 of file Particle.h.

◆ m_identifier

int m_identifier = -1
private

Identifier that can be used to identify whether the particle is unique or is a copy or representation of another.

For example a kaon and pion particles constructed from the same Track are representations of the same physical object in the detector and cannot be used in the reconstruction of the same decay chain

Definition at line 1084 of file Particle.h.

◆ m_jacobiMatrix

double m_jacobiMatrix[c_SizeMatrix] = {}
private

error matrix (1D representation)

Definition at line 1068 of file Particle.h.

◆ m_mass

double m_mass
private

particle (invariant) mass

Definition at line 1056 of file Particle.h.

◆ m_mdstIndex

unsigned m_mdstIndex
private

0-based index of MDST store array object

Definition at line 1073 of file Particle.h.

◆ m_momentumScale

double m_momentumScale = 1.0
private

effective momentum scale factor

Definition at line 1060 of file Particle.h.

◆ m_momentumScalingFactor

double m_momentumScalingFactor = 1.0
private

momentum scaling factor

Definition at line 1061 of file Particle.h.

◆ m_momentumSmearingFactor

double m_momentumSmearingFactor = 1.0
private

momentum smearing factor

Definition at line 1062 of file Particle.h.

◆ m_particleSource

EParticleSourceObject m_particleSource
private

(mdst) source of particle

Definition at line 1072 of file Particle.h.

◆ m_pdgCode

int m_pdgCode
private

PDG code.

Definition at line 1054 of file Particle.h.

◆ m_pdgCodeUsedForFit

int m_pdgCodeUsedForFit = 0
private

PDG code used for the track fit.

Definition at line 1055 of file Particle.h.

◆ m_properties

int m_properties
private

particle property

Definition at line 1074 of file Particle.h.

◆ m_pValue

double m_pValue
private

chi^2 probability of the fit.

Default is nan

Definition at line 1069 of file Particle.h.

◆ m_px

double m_px
private

momentum component x

Definition at line 1057 of file Particle.h.

◆ m_py

double m_py
private

momentum component y

Definition at line 1058 of file Particle.h.

◆ m_pz

double m_pz
private

momentum component z

Definition at line 1059 of file Particle.h.

◆ m_x

double m_x
private

position component x

Definition at line 1064 of file Particle.h.

◆ m_y

double m_y
private

position component y

Definition at line 1065 of file Particle.h.

◆ m_z

double m_z
private

position component z

Definition at line 1066 of file Particle.h.


The documentation for this class was generated from the following files: