Belle II Software  release-05-02-19
TrackFitResult.h
1 /**************************************************************************
2  * BASF2 (Belle Analysis Framework 2) *
3  * Copyright(C) 2013, 2014, 2015 - Belle II Collaboration *
4  * *
5  * Author: The Belle II Collaboration *
6  * Contributors: Martin Heck, Markus Prim, Tobias Schlüter *
7  * *
8  * This software is provided "as is" without any warranty. *
9  **************************************************************************/
10 #pragma once
11 
12 #include <framework/gearbox/Const.h>
13 #include <framework/datastore/RelationsObject.h>
14 #include <framework/dataobjects/Helix.h>
15 #include <framework/dataobjects/UncertainHelix.h>
16 #include <framework/geometry/BFieldManager.h>
17 
18 #include <TVector3.h>
19 #include <TMatrixDSym.h>
20 #include <TLorentzVector.h>
21 
22 #include <stdint.h>
23 #include <vector>
24 #include <cmath>
25 
26 namespace Belle2 {
31  class HitPatternCDC;
32  class HitPatternVXD;
33 
51  class TrackFitResult : public RelationsObject {
52  public:
55 
74  TrackFitResult(const TVector3& position, const TVector3& momentum,
75  const TMatrixDSym& covariance, const short int charge,
76  const Const::ParticleType& particleType, const float pValue,
77  const float bField,
78  const uint64_t hitPatternCDCInitializer,
79  const uint32_t hitPatternVXDInitializer,
80  const uint16_t NDF);
81 
93  TrackFitResult(const std::vector<float>& tau, const std::vector<float>& cov5,
94  const Const::ParticleType& particleType, const float pValue,
95  const uint64_t hitPatternCDCInitializer,
96  const uint32_t hitPatternVXDInitializer,
97  const uint16_t NDF
98  );
99 
101  TVector3 getPosition() const { return getHelix().getPerigee(); }
102 
108  TVector3 getMomentum() const
109  {
110  const double bField = BFieldManager::getField(getPosition()).Z() / Unit::T;
111  return getHelix().getMomentum(bField);
112  }
113 
117  TLorentzVector get4Momentum() const
118  {
119  return TLorentzVector(getMomentum(), getEnergy());
120  }
121 
125  double getEnergy() const
126  {
127  return std::sqrt(getMomentum().Mag2() + getParticleType().getMass() * getParticleType().getMass());
128  }
129 
132  double getTransverseMomentum() const
133  {
134  const double bField = BFieldManager::getField(getPosition()).Z() / Unit::T;
135  return getHelix().getTransverseMomentum(bField);
136  }
137 
139  TMatrixDSym getCovariance6() const
140  {
141  const double bField = BFieldManager::getField(getPosition()).Z() / Unit::T;
142  return getUncertainHelix().getCartesianCovariance(bField);
143  }
144 
147 
152  short getChargeSign() const { return (getOmega() > 0) - (getOmega() < 0); }
153 
155  double getPValue() const { return m_pValue; }
156 
158  int getNDF() const;
159 
161  double getChi2() const;
162 
163  //------------------------------------------------------------------------
164  // --- Getters for perigee helix parameters
165  //------------------------------------------------------------------------
170  double getD0() const { return m_tau[iD0]; }
171 
176  double getPhi0() const { return m_tau[iPhi0]; }
177 
179  double getPhi() const { return getPhi0(); }
180 
186  double getOmega() const { return m_tau[iOmega]; }
187 
192  double getZ0() const { return m_tau[iZ0]; }
193 
198  double getTanLambda() const { return m_tau[iTanLambda]; }
199 
201  double getCotTheta() const { return getTanLambda(); }
202 
207  std::vector<float> getTau() const { return std::vector<float>(m_tau, m_tau + c_NPars); }
208 
213  std::vector<float> getCov() const { return std::vector<float>(m_cov5, m_cov5 + c_NCovEntries); }
214 
219  TMatrixDSym getCovariance5() const;
220 
222  Helix getHelix() const
223  { return Helix(getD0(), getPhi0(), getOmega(), getZ0(), getTanLambda()); }
224 
226  UncertainHelix getUncertainHelix() const
227  {
228  return UncertainHelix(getD0(), getPhi0(), getOmega(), getZ0(),
230  }
231 
234 
237 
239  virtual std::string getInfoHTML() const override;
240 
242  private:
243 
244  //---------------------------------------------------------------------------------------------------------------------------
246  const unsigned int m_pdg;
247 
249  const Double32_t m_pValue;
250 
256  static const unsigned int c_NPars = 5;
257  static const unsigned int c_NCovEntries = 5 * 6 / 2;
258 
263  static const unsigned int iD0 = 0;
264  static const unsigned int iPhi0 = 1;
265  static const unsigned int iOmega = 2;
266  static const unsigned int iZ0 = 3;
267  static const unsigned int iTanLambda = 4;
268 
271  Double32_t m_tau[c_NPars];
272 
277  Double32_t m_cov5[c_NCovEntries];
278 
288  const uint32_t m_hitPatternVXDInitializer;
289 
291  static const uint16_t c_NDFFlag = 0xFFFF;
292 
294  uint16_t m_NDF;
295 
297  /* Version history:
298  ver 8: add NDF
299  ver 7: fixed sign errors in the translation of position and momentum covariances.
300  ver 6: use fixed size arrays instead of vectors (add schema evolution rule), use Double32_t.
301  ver 5: CDC Hit Pattern now a single variable (add schema evolution rule).
302  ver 4: added hit pattern
303  ver 3: back to vectors for the elements
304  ver 4: fixed-size arrays for the elements <------- incompatible with later version 4
305  ver 3: add data fields for hit patterns
306  ver 2: re-arranging of members
307  ver 1:
308  ver 2: <------- incompatible with later version 2
309  ver 1:
310  */
311  };
313 }
Belle2::TrackFitResult::TrackFitResult
TrackFitResult()
Constructor initializing everything to zero.
Definition: TrackFitResult.cc:24
Belle2::TrackFitResult::getNDF
int getNDF() const
Getter for number of degrees of freedom of the track fit.
Definition: TrackFitResult.cc:85
Belle2::TrackFitResult::get4Momentum
TLorentzVector get4Momentum() const
Getter for the 4Momentum at the closest approach of the track in the r/phi projection.
Definition: TrackFitResult.h:125
Belle2::TrackFitResult::iTanLambda
static const unsigned int iTanLambda
Index tan lambda.
Definition: TrackFitResult.h:275
Belle2::UncertainHelix::getCartesianCovariance
TMatrixDSym getCartesianCovariance(const double bZ_tesla=1.5) const
Getter for the position and momentum covariance matrix.
Definition: UncertainHelix.cc:115
Belle2::TrackFitResult::c_NDFFlag
static const uint16_t c_NDFFlag
backward compatibility initialisation for NDF
Definition: TrackFitResult.h:299
Belle2::BFieldManager::getField
static void getField(const double *pos, double *field)
return the magnetic field at a given position.
Definition: BFieldManager.h:110
Belle2::TrackFitResult::m_cov5
Double32_t m_cov5[c_NCovEntries]
The 15 = 5*(5+1)/2 covariance matrix elements.
Definition: TrackFitResult.h:285
Belle2::TrackFitResult::iD0
static const unsigned int iD0
Index for d0.
Definition: TrackFitResult.h:271
Belle2::TrackFitResult::getMomentum
TVector3 getMomentum() const
Getter for vector of momentum at closest approach of track in r/phi projection.
Definition: TrackFitResult.h:116
Belle2::TrackFitResult::getOmega
double getOmega() const
Getter for omega.
Definition: TrackFitResult.h:194
Belle2::TrackFitResult::getPValue
double getPValue() const
Getter for Chi2 Probability of the track fit.
Definition: TrackFitResult.h:163
Belle2::TrackFitResult::c_NCovEntries
static const unsigned int c_NCovEntries
Number of covariance entries.
Definition: TrackFitResult.h:265
Belle2::TrackFitResult::m_hitPatternVXDInitializer
const uint32_t m_hitPatternVXDInitializer
Member for initializing the information about hits in the VXD.
Definition: TrackFitResult.h:296
Belle2::TrackFitResult::m_pValue
const Double32_t m_pValue
Chi2 Probability of the fit.
Definition: TrackFitResult.h:257
Belle2::TrackFitResult::getTanLambda
double getTanLambda() const
Getter for tanLambda.
Definition: TrackFitResult.h:206
Belle2::TrackFitResult
Values of the result of a track fit with a given particle hypothesis.
Definition: TrackFitResult.h:59
Belle2::TrackFitResult::getPosition
TVector3 getPosition() const
Getter for vector of position at closest approach of track in r/phi projection.
Definition: TrackFitResult.h:109
Belle2::TrackFitResult::getHitPatternCDC
HitPatternCDC getHitPatternCDC() const
Getter for the hit pattern in the CDC;.
Definition: TrackFitResult.cc:120
Belle2::TrackFitResult::getZ0
double getZ0() const
Getter for z0.
Definition: TrackFitResult.h:200
Belle2::TrackFitResult::getCotTheta
double getCotTheta() const
Getter for tanLambda with CDF naming convention.
Definition: TrackFitResult.h:209
Belle2::Unit::T
static const double T
[tesla]
Definition: Unit.h:130
Belle2::TrackFitResult::getParticleType
Const::ParticleType getParticleType() const
Getter for ParticleType of the mass hypothesis of the track fit.
Definition: TrackFitResult.h:154
Belle2::TrackFitResult::m_NDF
uint16_t m_NDF
Memeber for number of degrees of freedom.
Definition: TrackFitResult.h:302
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:19
Belle2::TrackFitResult::iPhi0
static const unsigned int iPhi0
Index for phi0.
Definition: TrackFitResult.h:272
Belle2::TrackFitResult::getCov
std::vector< float > getCov() const
Getter for all covariance matrix elements of perigee parameters.
Definition: TrackFitResult.h:221
Belle2::TrackFitResult::m_hitPatternCDCInitializer
uint64_t m_hitPatternCDCInitializer
Member for initializing the information about hits in the CDC.
Definition: TrackFitResult.h:291
Belle2::TrackFitResult::getTransverseMomentum
double getTransverseMomentum() const
Getter for the absolute value of the transverse momentum at the perigee.
Definition: TrackFitResult.h:140
Belle2::TrackFitResult::getPhi0
double getPhi0() const
Getter for phi0.
Definition: TrackFitResult.h:184
Belle2::TrackFitResult::getChi2
double getChi2() const
Get chi2 given NDF and p-value.
Definition: TrackFitResult.cc:93
Belle2::TrackFitResult::getHitPatternVXD
HitPatternVXD getHitPatternVXD() const
Getter for the hit pattern in the VXD;.
Definition: TrackFitResult.cc:125
Belle2::TrackFitResult::m_pdg
const unsigned int m_pdg
PDG Code for hypothesis with which the corresponding fit was performed.
Definition: TrackFitResult.h:254
Belle2::TrackFitResult::getUncertainHelix
UncertainHelix getUncertainHelix() const
Conversion to framework Uncertain Helix (i.e., with covariance).
Definition: TrackFitResult.h:234
Belle2::TrackFitResult::getCovariance5
TMatrixDSym getCovariance5() const
Getter for covariance matrix of perigee parameters in matrix form.
Definition: TrackFitResult.cc:106
Belle2::Const::ParticleType
The ParticleType class for identifying different particle types.
Definition: Const.h:284
Belle2::TrackFitResult::iOmega
static const unsigned int iOmega
Index for omega.
Definition: TrackFitResult.h:273
Belle2::RelationsObject
RelationsInterface< TObject > RelationsObject
Provides interface for getting/adding relations to objects in StoreArrays.
Definition: RelationsObject.h:443
Belle2::TrackFitResult::getInfoHTML
virtual std::string getInfoHTML() const override
Return a short summary of this object's contents in HTML format.
Definition: TrackFitResult.cc:129
Belle2::TrackFitResult::ClassDefOverride
ClassDefOverride(TrackFitResult, 8)
Values of the result of a track fit with a given particle hypothesis.
Belle2::TrackFitResult::getHelix
Helix getHelix() const
Conversion to framework Helix (without covariance).
Definition: TrackFitResult.h:230
Belle2::TrackFitResult::c_NPars
static const unsigned int c_NPars
Number of helix parameters.
Definition: TrackFitResult.h:264
Belle2::TrackFitResult::getTau
std::vector< float > getTau() const
Getter for all perigee parameters.
Definition: TrackFitResult.h:215
Belle2::TrackFitResult::getCovariance6
TMatrixDSym getCovariance6() const
Position and Momentum Covariance Matrix.
Definition: TrackFitResult.h:147
Belle2::TrackFitResult::getEnergy
double getEnergy() const
Getter for the Energy at the closest approach of the track in the r/phi projection.
Definition: TrackFitResult.h:133
Belle2::HitPatternVXD
Hit pattern of the VXD within a track.
Definition: HitPatternVXD.h:44
Belle2::TrackFitResult::getChargeSign
short getChargeSign() const
Return track charge (1 or -1).
Definition: TrackFitResult.h:160
Belle2::HitPatternCDC
Hit pattern of CDC hits within a track.
Definition: HitPatternCDC.h:45
Belle2::TrackFitResult::getPhi
double getPhi() const
Getter for phi0 with CDF naming convention.
Definition: TrackFitResult.h:187
Belle2::TrackFitResult::m_tau
Double32_t m_tau[c_NPars]
perigee helix parameters; tau = d0, phi0, omega, z0, tanLambda.
Definition: TrackFitResult.h:279
Belle2::TrackFitResult::getD0
double getD0() const
Getter for d0.
Definition: TrackFitResult.h:178
Belle2::TrackFitResult::iZ0
static const unsigned int iZ0
Index for z0.
Definition: TrackFitResult.h:274