ATLAS Offline Software
Loading...
Searching...
No Matches
xAOD::NeutralParticle_v1 Class Reference

Class describing a NeutralParticle. More...

#include <NeutralParticle_v1.h>

Inheritance diagram for xAOD::NeutralParticle_v1:
Collaboration diagram for xAOD::NeutralParticle_v1:

Public Member Functions

 NeutralParticle_v1 ()
 Default constructor.
 ~NeutralParticle_v1 ()
 Destructor.
 NeutralParticle_v1 (const NeutralParticle_v1 &o)
 Copy ctor. This involves copying the entire Auxilary store, and is a slow operation which should be used sparingly.
NeutralParticle_v1operator= (const NeutralParticle_v1 &tp)
 Assignment operator. This can involve creating and copying an Auxilary store, and so should be used sparingly.

xAOD::IParticle functions

The transverse momentum ( \(p_T\)) of the particle.

typedef IParticle::FourMom_t FourMom_t
 Definition of the 4-momentum type.
typedef ROOT::Math::LorentzVector< ROOT::Math::PxPyPzM4D< double > > GenVecFourMom_t
 Base 4 Momentum type for TrackParticle.
virtual double pt () const
 The transverse momentum ( \(p_T\)) of the particle.
virtual double eta () const
 The pseudorapidity ( \(\eta\)) of the particle.
virtual double phi () const
 The azimuthal angle ( \(\phi\)) of the particle.
virtual double m () const
 The invariant mass of the particle..
virtual double e () const
 The total energy of the particle.
virtual double rapidity () const
 The true rapidity (y) of the particle.
virtual FourMom_t p4 () const
 The full 4-momentum of the particle.
GenVecFourMom_t genvecP4 () const
 The full 4-momentum of the particle : GenVector form.
virtual Type::ObjectType type () const
 The type of the object as a simple enumeration.

Defining parameters functions

The 'defining parameters' are key to the concept of a NeutralParticle, and give the values for the IParticle interface ( pt(), phi(), eta() etc.).

They use the Trk::NeutralPerigee coordinate system, and are defined as: \(( d_0, z_0, \phi, \theta, q/p )\). The parameters are expressed with respect to an origin (returned by vx(), vy() and vy() ), currently intended to be the 'beamspot'. This origin is expected to be the same for all track particles in a collection (and this may be be enforced).

CxxUtils::CachedValue< Trk::NeutralPerigeem_perigeeParameters
 Cached NeutralPerigee, built from this object.
float d0 () const
 Returns the \(d_0\) parameter.
float z0 () const
 Returns the \(z_0\) parameter.
float phi0 () const
 Returns the \(\phi\) parameter, which has range \(-\pi\) to \(+\pi\).
float theta () const
 Returns the \(\theta\) parameter, which has range 0 to \(\pi\).
float oneOverP () const
 Returns the \(q/p\) parameter.
const DefiningParameters_t definingParameters () const
 Returns a SVector of the Perigee track parameters.
const ParametersCovMatrix_t definingParametersCovMatrix () const
 Returns the 5x5 symmetric matrix containing the defining parameters covariance matrix.
const std::vector< float > & definingParametersCovMatrixVec () const
 Returns the vector of the covariance values - 15 elements.
void setDefiningParameters (float d0, float z0, float phi0, float theta, float qOverP)
 Set the defining parameters.
void setDefiningParametersCovMatrix (const ParametersCovMatrix_t &cov)
 Set the defining parameters covariance matrix.
void setDefiningParametersCovMatrixVec (const std::vector< float > &cov)
 Set the defining parameters covariance matrix using a length 15 vector.
float vx () const
 The x origin for the parameters.
float vy () const
 The y origin for the parameters.
float vz () const
 The z origin for the parameters.
void setParametersOrigin (float x, float y, float z)
 Set the origin for the parameters.
const Trk::NeutralPerigeeperigeeParameters () const
 Returns the Trk::NeutralPerigee track parameters.
void resetCache ()
 Reset the internal cache of the object.

Functions for getting and setting user properties

template<class T>
XAOD_AUXDATA_DEPRECATED T & auxdata (const std::string &name, const std::string &clsname="")
 Fetch an aux data variable, as a non-const reference.
template<class T>
XAOD_AUXDATA_DEPRECATED const T & auxdata (const std::string &name, const std::string &clsname="") const
 Fetch an aux data variable, as a const reference.
template<class T>
XAOD_AUXDATA_DEPRECATED bool isAvailable (const std::string &name, const std::string &clsname="") const
 Check if a user property is available for reading or not.
template<class T>
XAOD_AUXDATA_DEPRECATED bool isAvailableWritable (const std::string &name, const std::string &clsname="") const
 Check if a user property is available for writing or not.

Detailed Description

Class describing a NeutralParticle.

Author
Andreas Salzburger Andre.nosp@m.as.S.nosp@m.alzbu.nosp@m.rger.nosp@m.@cern.nosp@m..ch

Definition at line 40 of file NeutralParticle_v1.h.

Member Typedef Documentation

◆ FourMom_t

Definition of the 4-momentum type.

Definition at line 69 of file NeutralParticle_v1.h.

◆ GenVecFourMom_t

typedef ROOT::Math::LorentzVector<ROOT::Math::PxPyPzM4D<double> > xAOD::NeutralParticle_v1::GenVecFourMom_t

Base 4 Momentum type for TrackParticle.

Definition at line 75 of file NeutralParticle_v1.h.

Constructor & Destructor Documentation

◆ NeutralParticle_v1() [1/2]

xAOD::NeutralParticle_v1::NeutralParticle_v1 ( )

Default constructor.

Definition at line 20 of file NeutralParticle_v1.cxx.

21 : IParticle(){
22 }
IParticle()=default

◆ ~NeutralParticle_v1()

xAOD::NeutralParticle_v1::~NeutralParticle_v1 ( )

Destructor.

Definition at line 40 of file NeutralParticle_v1.cxx.

40 {
41 }

◆ NeutralParticle_v1() [2/2]

xAOD::NeutralParticle_v1::NeutralParticle_v1 ( const NeutralParticle_v1 & o)

Copy ctor. This involves copying the entire Auxilary store, and is a slow operation which should be used sparingly.

Definition at line 24 of file NeutralParticle_v1.cxx.

25 : IParticle( tp ) {
26 makePrivateStore( tp );
27 }
void makePrivateStore()
Create a new (empty) private store for this object.

Member Function Documentation

◆ auxdata() [1/2]

template<class T>
XAOD_AUXDATA_DEPRECATED T & xAOD::IParticle::auxdata ( const std::string & name,
const std::string & clsname = "" )
inlineinherited

Fetch an aux data variable, as a non-const reference.

This function provides an easy way for users to decorate objects with auxiliary data.

Take note that this function is slow. Should not be used inside time-critical code.

Parameters
nameName of the aux variable
clsnameThe name of the associated class. May be blank
Returns
A modifiable reference to the decoration

Definition at line 98 of file Event/xAOD/xAODBase/xAODBase/IParticle.h.

99 {
100
101 return SG::Accessor< T >(name, clsname)(*this);
102 }
SG::Accessor< T, ALLOC > Accessor
Definition AuxElement.h:573

◆ auxdata() [2/2]

template<class T>
XAOD_AUXDATA_DEPRECATED const T & xAOD::IParticle::auxdata ( const std::string & name,
const std::string & clsname = "" ) const
inlineinherited

Fetch an aux data variable, as a const reference.

This function provides an easy way for users to retrieve auxiliary decorations from an object.

Take note that this function is slow. Should not be used inside time-critical code.

Parameters
nameName of the aux variable
clsnameThe name of the associated class. May be blank
Returns
A constant reference to the decoration

Definition at line 118 of file Event/xAOD/xAODBase/xAODBase/IParticle.h.

119 {
120
121 return SG::ConstAccessor< T >( name, clsname )( *this );
122 }
SG::ConstAccessor< T, ALLOC > ConstAccessor
Definition AuxElement.h:570

◆ d0()

float xAOD::NeutralParticle_v1::d0 ( ) const

Returns the \(d_0\) parameter.

◆ definingParameters()

const DefiningParameters_t xAOD::NeutralParticle_v1::definingParameters ( ) const

Returns a SVector of the Perigee track parameters.

i.e. a vector of \(\left(\begin{array}{c}d_0\\z_0\\\phi_0\\\theta\\q/p\end{array}\right)\)

Definition at line 104 of file NeutralParticle_v1.cxx.

104 {
105 DefiningParameters_t tmp;
106 tmp << d0(),z0(),phi0(),theta(),oneOverP();
107 return tmp;
108 }
float d0() const
Returns the parameter.
float oneOverP() const
Returns the parameter.
float theta() const
Returns the parameter, which has range 0 to .
float phi0() const
Returns the parameter, which has range to .
float z0() const
Returns the parameter.

◆ definingParametersCovMatrix()

const xAOD::ParametersCovMatrix_t xAOD::NeutralParticle_v1::definingParametersCovMatrix ( ) const

Returns the 5x5 symmetric matrix containing the defining parameters covariance matrix.

Definition at line 149 of file NeutralParticle_v1.cxx.

149 {
150 static const Accessor< std::vector<float> > acc( "definingParametersCovMatrix" );
151 std::vector<float> v = acc(*this);
152 std::vector<float>::const_iterator it = v.begin();
153 xAOD::ParametersCovMatrix_t cov;
154 cov.setZero();
155 for (size_t irow = 0; irow<5; ++irow){
156 for (size_t icol =0; icol<=irow; ++icol){
157 cov.fillSymmetric(icol,irow, *it++);
158 }
159 }
160 return cov;
161 }
static const SG::AuxElement::Accessor< ElementLink< IParticleContainer > > acc("originalObjectLink")
Object used for setting/getting the dynamic decoration in question.

◆ definingParametersCovMatrixVec()

const std::vector< float > & xAOD::NeutralParticle_v1::definingParametersCovMatrixVec ( ) const

Returns the vector of the covariance values - 15 elements.

Definition at line 163 of file NeutralParticle_v1.cxx.

163 {
164 static const Accessor< std::vector<float> > acc( "definingParametersCovMatrix" );
165 return acc(*this);
166 }

◆ e()

virtual double xAOD::NeutralParticle_v1::e ( ) const
virtual

The total energy of the particle.

Implements xAOD::IParticle.

◆ eta()

double xAOD::NeutralParticle_v1::eta ( ) const
virtual

The pseudorapidity ( \(\eta\)) of the particle.

Implements xAOD::IParticle.

Definition at line 47 of file NeutralParticle_v1.cxx.

47 {
48 return genvecP4().Eta();
49 }
GenVecFourMom_t genvecP4() const
The full 4-momentum of the particle : GenVector form.

◆ genvecP4()

NeutralParticle_v1::GenVecFourMom_t xAOD::NeutralParticle_v1::genvecP4 ( ) const

The full 4-momentum of the particle : GenVector form.

Definition at line 65 of file NeutralParticle_v1.cxx.

65 {
66 using namespace std;
67 float p = 1/fabs(oneOverP());
68 float thetaT = theta();
69 float phiT = phi();
70 float sinTheta= sin(thetaT);
71 float px = p*sinTheta*cos(phiT);
72 float py = p*sinTheta*sin(phiT);
73 float pz = p*cos(thetaT);
74 return GenVecFourMom_t(px, py, pz, m());
75 }
virtual double m() const
The invariant mass of the particle..
ROOT::Math::LorentzVector< ROOT::Math::PxPyPzM4D< double > > GenVecFourMom_t
Base 4 Momentum type for TrackParticle.

◆ isAvailable()

template<class T>
XAOD_AUXDATA_DEPRECATED bool xAOD::IParticle::isAvailable ( const std::string & name,
const std::string & clsname = "" ) const
inlineinherited

Check if a user property is available for reading or not.

This function should be used to check if a user property which may or may not exist, is set on the object.

Parameters
nameName of the auxiliary variable
clsnameThe name of the associated class. May be blank
Returns
Whether the decoration exists or not

Definition at line 135 of file Event/xAOD/xAODBase/xAODBase/IParticle.h.

136 {
137
138 return SG::ConstAccessor< T >(name, clsname).isAvailable(*this);
139 }
bool isAvailable(const ELT &e) const
Test to see if this variable exists in the store.

◆ isAvailableWritable()

template<class T>
XAOD_AUXDATA_DEPRECATED bool xAOD::IParticle::isAvailableWritable ( const std::string & name,
const std::string & clsname = "" ) const
inlineinherited

Check if a user property is available for writing or not.

This function can be used to check whether it will be possible to set a user property on the object.

Parameters
nameName of the auxiliary variable
clsnameThe name of the associated class. May be blank
Returns
Whether the decoration is possible to set

Definition at line 152 of file Event/xAOD/xAODBase/xAODBase/IParticle.h.

153 {
154
155 return SG::Accessor< T >(name, clsname).isAvailableWritable(*this);
156 }
bool isAvailableWritable(ELT &e) const
Test to see if this variable exists in the store and is writable.

◆ m()

double xAOD::NeutralParticle_v1::m ( ) const
virtual

The invariant mass of the particle..

Todo
Get value from somewhere. Also, the TrackParticle took the Pion mass - do we really want to do this? We have ParticleHypo?

Implements xAOD::IParticle.

Definition at line 54 of file NeutralParticle_v1.cxx.

54 {
56 }
constexpr double chargedPionMassInMeV
the mass of the charged pion (in MeV)

◆ oneOverP()

float xAOD::NeutralParticle_v1::oneOverP ( ) const

Returns the \(q/p\) parameter.

◆ operator=()

NeutralParticle_v1 & xAOD::NeutralParticle_v1::operator= ( const NeutralParticle_v1 & tp)

Assignment operator. This can involve creating and copying an Auxilary store, and so should be used sparingly.

Definition at line 29 of file NeutralParticle_v1.cxx.

29 {
30 if(this == &tp) return *this;
31
32 if(!hasStore() ) makePrivateStore();
33 this->IParticle::operator=( tp );
34#ifndef XAOD_ANALYSIS
35 m_perigeeParameters.reset();
36#endif // not XAOD_ANALYSIS
37 return *this;
38 }
IParticle & operator=(const IParticle &)=default
CxxUtils::CachedValue< Trk::NeutralPerigee > m_perigeeParameters
Cached NeutralPerigee, built from this object.
bool hasStore() const
Return true if this object has an associated store.

◆ p4()

NeutralParticle_v1::FourMom_t xAOD::NeutralParticle_v1::p4 ( ) const
virtual

The full 4-momentum of the particle.

Implements xAOD::IParticle.

Definition at line 77 of file NeutralParticle_v1.cxx.

77 {
78 using namespace std;
80 float p = 1/fabs(oneOverP());
81 float thetaT = theta();
82 float phiT = phi();
83 float sinTheta= sin(thetaT);
84 float px = p*sinTheta*cos(phiT);
85 float py = p*sinTheta*sin(phiT);
86 float pz = p*cos(thetaT);
87 float e = pow (m(),2) +
88 pow( px,2) + pow( py,2) + pow( pz,2);
89 p4.SetPxPyPzE( px, py, pz, sqrt(e) );
90 return p4;
91 }
IParticle::FourMom_t FourMom_t
Definition of the 4-momentum type.
virtual FourMom_t p4() const
The full 4-momentum of the particle.
virtual double e() const
The total energy of the particle.
constexpr int pow(int x)
Definition conifer.h:27

◆ perigeeParameters()

const Trk::NeutralPerigee & xAOD::NeutralParticle_v1::perigeeParameters ( ) const

Returns the Trk::NeutralPerigee track parameters.

These are defined as: \(\left(\begin{array}{c}d_0\\z_0\\\phi_0\\\theta\\1/p\\\end{array}\right)\)

Note
This is only available in Athena.

Definition at line 189 of file NeutralParticle_v1.cxx.

189 {
190
191 // Require the cache to be valid and check if the cached pointer has been set
192 if(m_perigeeParameters.isValid()){
193 return *(m_perigeeParameters.ptr());
194 }
195 static const Accessor< float > acc1( "d0" );
196 static const Accessor< float > acc2( "z0" );
197 static const Accessor< float > acc3( "phi" );
198 static const Accessor< float > acc4( "theta" );
199 static const Accessor< float > acc5( "oneOverP" );
200 static const Accessor< std::vector<float> > acc6( "definingParametersCovMatrix" );
201 ParametersCovMatrix_t cov;
202 cov.setZero();
203 auto it= acc6(*this).begin();
204 for (size_t irow = 0; irow<5; ++irow){
205 for (size_t icol =0; icol<=irow; ++icol){
206 cov.fillSymmetric(irow,icol,*it++) ;
207 }
208 }
209 Trk::NeutralPerigee tmpPerigeeParameters(acc1(*this),acc2(*this),acc3(*this),acc4(*this),acc5(*this),
210 Trk::PerigeeSurface(Amg::Vector3D(vx(),vy(),vz())),std::move(cov));
211 m_perigeeParameters.set(tmpPerigeeParameters);
212 return *(m_perigeeParameters.ptr());
213 }
float vz() const
The z origin for the parameters.
float vx() const
The x origin for the parameters.
float vy() const
The y origin for the parameters.
Eigen::Matrix< double, 3, 1 > Vector3D
ParametersT< NeutralParametersDim, Neutral, PerigeeSurface > NeutralPerigee

◆ phi()

virtual double xAOD::NeutralParticle_v1::phi ( ) const
virtual

The azimuthal angle ( \(\phi\)) of the particle.

Implements xAOD::IParticle.

◆ phi0()

float xAOD::NeutralParticle_v1::phi0 ( ) const

Returns the \(\phi\) parameter, which has range \(-\pi\) to \(+\pi\).

◆ pt()

double xAOD::NeutralParticle_v1::pt ( ) const
virtual

The transverse momentum ( \(p_T\)) of the particle.

Implements xAOD::IParticle.

Definition at line 43 of file NeutralParticle_v1.cxx.

43 {
44 return genvecP4().Pt();
45 }

◆ rapidity()

double xAOD::NeutralParticle_v1::rapidity ( ) const
virtual

The true rapidity (y) of the particle.

Implements xAOD::IParticle.

Definition at line 61 of file NeutralParticle_v1.cxx.

61 {
62 return genvecP4().Rapidity();
63 }

◆ resetCache()

void xAOD::NeutralParticle_v1::resetCache ( )

Reset the internal cache of the object.

Definition at line 216 of file NeutralParticle_v1.cxx.

216 {
217#ifndef XAOD_ANALYSIS
218 m_perigeeParameters.reset();
219#endif // not XAOD_ANALYSIS
220 }

◆ setDefiningParameters()

void xAOD::NeutralParticle_v1::setDefiningParameters ( float d0,
float z0,
float phi0,
float theta,
float qOverP )

Set the defining parameters.

Definition at line 110 of file NeutralParticle_v1.cxx.

110 {
111#ifndef XAOD_ANALYSIS
112 if(m_perigeeParameters.isValid()) {
113 m_perigeeParameters.reset();
114 }
115#endif // not XAOD_ANALYSIS
116 static const Accessor< float > acc1( "d0" );
117 acc1( *this ) = d0;
118
119 static const Accessor< float > acc2( "z0" );
120 acc2( *this ) = z0;
121
122 static const Accessor< float > acc3( "phi" );
123 acc3( *this ) = phi0;
124
125 static const Accessor< float > acc4( "theta" );
126 acc4( *this ) = theta;
127
128 static const Accessor< float > acc5( "oneOverP" );
129 acc5( *this ) = oneOverP;
130
131 return;
132 }

◆ setDefiningParametersCovMatrix()

void xAOD::NeutralParticle_v1::setDefiningParametersCovMatrix ( const ParametersCovMatrix_t & cov)

Set the defining parameters covariance matrix.

Definition at line 134 of file NeutralParticle_v1.cxx.

134 {
135#ifndef XAOD_ANALYSIS
136 if(m_perigeeParameters.isValid()) {
137 m_perigeeParameters.reset();
138 }
139#endif // not XAOD_ANALYSIS
140
141 static const Accessor< std::vector<float> > acc( "definingParametersCovMatrix" );
142 std::vector<float>& v = acc(*this);
143 v.reserve(15);
144 for (size_t irow = 0; irow<5; ++irow)
145 for (size_t icol =0; icol<=irow; ++icol)
146 v.push_back(cov(icol,irow));
147 }

◆ setDefiningParametersCovMatrixVec()

void xAOD::NeutralParticle_v1::setDefiningParametersCovMatrixVec ( const std::vector< float > & cov)

Set the defining parameters covariance matrix using a length 15 vector.

Definition at line 168 of file NeutralParticle_v1.cxx.

168 {
169 static const Accessor< std::vector<float> > acc( "definingParametersCovMatrix" );
170 acc(*this)=cov;
171 }

◆ setParametersOrigin()

void xAOD::NeutralParticle_v1::setParametersOrigin ( float x,
float y,
float z )

Set the origin for the parameters.

Definition at line 177 of file NeutralParticle_v1.cxx.

177 {
178 static const Accessor< float > acc1( "vx" );
179 acc1( *this ) = x;
180
181 static const Accessor< float > acc2( "vy" );
182 acc2( *this ) = y;
183
184 static const Accessor< float > acc3( "vz" );
185 acc3( *this ) = z;
186 }
#define y
#define x
#define z

◆ theta()

float xAOD::NeutralParticle_v1::theta ( ) const

Returns the \(\theta\) parameter, which has range 0 to \(\pi\).

◆ type()

Type::ObjectType xAOD::NeutralParticle_v1::type ( ) const
virtual

The type of the object as a simple enumeration.

Implements xAOD::IParticle.

Definition at line 93 of file NeutralParticle_v1.cxx.

93 {
94 return Type::NeutralParticle;
95 }

◆ vx()

float xAOD::NeutralParticle_v1::vx ( ) const

The x origin for the parameters.

◆ vy()

float xAOD::NeutralParticle_v1::vy ( ) const

The y origin for the parameters.

◆ vz()

float xAOD::NeutralParticle_v1::vz ( ) const

The z origin for the parameters.

◆ z0()

float xAOD::NeutralParticle_v1::z0 ( ) const

Returns the \(z_0\) parameter.

Member Data Documentation

◆ m_perigeeParameters

CxxUtils::CachedValue<Trk::NeutralPerigee> xAOD::NeutralParticle_v1::m_perigeeParameters
private

Cached NeutralPerigee, built from this object.

Note
This is only available in Athena.

Definition at line 144 of file NeutralParticle_v1.h.


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