Belle II Software  light-2212-foldex
TrackFitResult.h
1 /**************************************************************************
2  * basf2 (Belle II Analysis Software Framework) *
3  * Author: The Belle II Collaboration *
4  * *
5  * See git log for contributors and copyright holders. *
6  * This file is licensed under LGPL-3.0, see LICENSE.md. *
7  **************************************************************************/
8 #pragma once
9 
10 #include <framework/gearbox/Const.h>
11 #include <framework/datastore/RelationsObject.h>
12 #include <framework/dataobjects/Helix.h>
13 #include <framework/dataobjects/UncertainHelix.h>
14 #include <framework/geometry/BFieldManager.h>
15 
16 #include <TMatrixDSym.h>
17 #include <Math/Vector3D.h>
18 #include <Math/Vector4D.h>
19 
20 #include <stdint.h>
21 #include <vector>
22 #include <cmath>
23 
24 namespace Belle2 {
29  class HitPatternCDC;
30  class HitPatternVXD;
31 
50  public:
53 
72  TrackFitResult(const ROOT::Math::XYZVector& position, const ROOT::Math::XYZVector& momentum,
73  const TMatrixDSym& covariance, const short int charge,
74  const Const::ParticleType& particleType, const float pValue,
75  const float bField,
76  const uint64_t hitPatternCDCInitializer,
77  const uint32_t hitPatternVXDInitializer,
78  const uint16_t NDF);
79 
91  TrackFitResult(const std::vector<float>& tau, const std::vector<float>& cov5,
92  const Const::ParticleType& particleType, const float pValue,
93  const uint64_t hitPatternCDCInitializer,
94  const uint32_t hitPatternVXDInitializer,
95  const uint16_t NDF
96  );
101  void updateTrackFitResult(const TrackFitResult& input);
102 
104  void mask() {m_pValue = NAN;}
105 
107  ROOT::Math::XYZVector getPosition() const { return getHelix().getPerigee(); }
108 
114  ROOT::Math::XYZVector getMomentum() const
115  {
116  const double bField = BFieldManager::getField(getPosition()).Z() / Unit::T;
117  return getHelix().getMomentum(bField);
118  }
119 
123  ROOT::Math::PxPyPzEVector get4Momentum() const
124  {
125  const ROOT::Math::XYZVector momentum = getMomentum();
126  return ROOT::Math::PxPyPzEVector(momentum.X(), momentum.Y(), momentum.Z(), getEnergy());
127  }
128 
132  double getEnergy() const
133  {
134  return std::sqrt(getMomentum().Mag2() + getParticleType().getMass() * getParticleType().getMass());
135  }
136 
139  double getTransverseMomentum() const
140  {
141  const double bField = BFieldManager::getField(getPosition()).Z() / Unit::T;
142  return getHelix().getTransverseMomentum(bField);
143  }
144 
146  TMatrixDSym getCovariance6() const
147  {
148  const double bField = BFieldManager::getField(getPosition()).Z() / Unit::T;
149  return getUncertainHelix().getCartesianCovariance(bField);
150  }
151 
154 
159  short getChargeSign() const { return (getOmega() > 0) - (getOmega() < 0); }
160 
162  double getPValue() const { return m_pValue; }
163 
165  int getNDF() const;
166 
168  double getChi2() const;
169 
170  //------------------------------------------------------------------------
171  // --- Getters for perigee helix parameters
172  //------------------------------------------------------------------------
177  double getD0() const { return m_tau[iD0]; }
178 
183  double getPhi0() const { return m_tau[iPhi0]; }
184 
186  double getPhi() const { return getPhi0(); }
187 
193  double getOmega() const { return m_tau[iOmega]; }
194 
199  double getZ0() const { return m_tau[iZ0]; }
200 
205  double getTanLambda() const { return m_tau[iTanLambda]; }
206 
208  double getCotTheta() const { return getTanLambda(); }
209 
214  std::vector<float> getTau() const { return std::vector<float>(m_tau, m_tau + c_NPars); }
215 
220  std::vector<float> getCov() const { return std::vector<float>(m_cov5, m_cov5 + c_NCovEntries); }
221 
226  TMatrixDSym getCovariance5() const;
227 
229  Helix getHelix() const
230  { return Helix(getD0(), getPhi0(), getOmega(), getZ0(), getTanLambda()); }
231 
233  Helix getHelix(float momentumScale) const
234  { return Helix(getD0(), getPhi0(), getOmega() / momentumScale, getZ0(), getTanLambda()); }
235 
238  {
239  return UncertainHelix(getD0(), getPhi0(), getOmega(), getZ0(),
241  }
242 
245 
248 
250  virtual std::string getInfoHTML() const override;
251 
253  private:
254 
255  //---------------------------------------------------------------------------------------------------------------------------
257  unsigned int m_pdg;
258 
260  Double32_t m_pValue;
261 
267  static const unsigned int c_NPars = 5;
268  static const unsigned int c_NCovEntries = 5 * 6 / 2;
274  static const unsigned int iD0 = 0;
275  static const unsigned int iPhi0 = 1;
276  static const unsigned int iOmega = 2;
277  static const unsigned int iZ0 = 3;
278  static const unsigned int iTanLambda = 4;
282  Double32_t m_tau[c_NPars];
283 
288  Double32_t m_cov5[c_NCovEntries];
289 
300 
302  static const uint16_t c_NDFFlag = 0xFFFF;
303 
305  uint16_t m_NDF;
306 
308  /* Version history:
309  ver 9: change m_pValue, m_pdg, m_hitPatternVXDInitializer to a non-const value
310  ver 8: add NDF
311  ver 7: fixed sign errors in the translation of position and momentum covariances.
312  ver 6: use fixed size arrays instead of vectors (add schema evolution rule), use Double32_t.
313  ver 5: CDC Hit Pattern now a single variable (add schema evolution rule).
314  ver 4: added hit pattern
315  ver 3: back to vectors for the elements
316  ver 4: fixed-size arrays for the elements <------- incompatible with later version 4
317  ver 3: add data fields for hit patterns
318  ver 2: re-arranging of members
319  ver 1:
320  ver 2: <------- incompatible with later version 2
321  ver 1:
322  */
323  friend class PostMergeUpdaterModule;
324  };
326 }
The ParticleType class for identifying different particle types.
Definition: Const.h:399
This class represents an ideal helix in perigee parameterization.
Definition: Helix.h:59
ROOT::Math::XYZVector getMomentum(const double bZ) const
Getter for vector of momentum at the perigee position.
Definition: Helix.cc:89
ROOT::Math::XYZVector getPerigee() const
Getter for the perigee position.
Definition: Helix.cc:69
double getTransverseMomentum(const double bZ) const
Getter for the absolute value of the transverse momentum at the perigee.
Definition: Helix.cc:99
Hit pattern of CDC hits within a track.
Definition: HitPatternCDC.h:35
Hit pattern of the VXD within a track.
Definition: HitPatternVXD.h:37
If the content of two DataStores are merged using the 'MergeDataStoreModule', several kinematic prope...
Defines interface for accessing relations of objects in StoreArray.
Values of the result of a track fit with a given particle hypothesis.
std::vector< float > getTau() const
Getter for all perigee parameters.
static const unsigned int c_NPars
Number of helix parameters.
Helix getHelix() const
Conversion to framework Helix (without covariance).
TMatrixDSym getCovariance5() const
Getter for covariance matrix of perigee parameters in matrix form.
uint32_t m_hitPatternVXDInitializer
Member for initializing the information about hits in the VXD.
static const unsigned int iZ0
Index for z0.
static const unsigned int c_NCovEntries
Number of covariance entries.
Helix getHelix(float momentumScale) const
Conversion to framework Helix with momentum scaling (without covariance).
double getPhi() const
Getter for phi0 with CDF naming convention.
unsigned int m_pdg
PDG Code for hypothesis with which the corresponding fit was performed.
Double32_t m_tau[c_NPars]
perigee helix parameters; tau = d0, phi0, omega, z0, tanLambda.
double getChi2() const
Get chi2 given NDF and p-value.
short getChargeSign() const
Return track charge (1 or -1).
double getPValue() const
Getter for Chi2 Probability of the track fit.
Double32_t m_pValue
Chi2 Probability of the fit.
double getEnergy() const
Getter for the Energy at the closest approach of the track in the r/phi projection.
TMatrixDSym getCovariance6() const
Position and Momentum Covariance Matrix.
double getCotTheta() const
Getter for tanLambda with CDF naming convention.
virtual std::string getInfoHTML() const override
Return a short summary of this object's contents in HTML format.
ROOT::Math::PxPyPzEVector get4Momentum() const
Getter for the 4Momentum at the closest approach of the track in the r/phi projection.
double getOmega() const
Getter for omega.
void updateTrackFitResult(const TrackFitResult &input)
update the TrackFitResults
int getNDF() const
Getter for number of degrees of freedom of the track fit.
double getD0() const
Getter for d0.
double getTransverseMomentum() const
Getter for the absolute value of the transverse momentum at the perigee.
double getTanLambda() const
Getter for tanLambda.
static const unsigned int iTanLambda
Index tan lambda.
TrackFitResult()
Constructor initializing everything to zero.
std::vector< float > getCov() const
Getter for all covariance matrix elements of perigee parameters.
static const unsigned int iD0
Index for d0.
uint64_t m_hitPatternCDCInitializer
Member for initializing the information about hits in the CDC.
void mask()
mask this TrackFitResults by setting the pValue to nan
static const unsigned int iOmega
Index for omega.
Const::ParticleType getParticleType() const
Getter for ParticleType of the mass hypothesis of the track fit.
static const uint16_t c_NDFFlag
backward compatibility initialisation for NDF
double getZ0() const
Getter for z0.
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.
uint16_t m_NDF
Memeber for number of degrees of freedom.
Double32_t m_cov5[c_NCovEntries]
The 15 = 5*(5+1)/2 covariance matrix elements.
static const unsigned int iPhi0
Index for phi0.
ClassDefOverride(TrackFitResult, 9)
Values of the result of a track fit with a given particle hypothesis.
HitPatternCDC getHitPatternCDC() const
Getter for the hit pattern in the CDC;.
UncertainHelix getUncertainHelix() const
Conversion to framework Uncertain Helix (i.e., with covariance).
double getPhi0() const
Getter for phi0.
HitPatternVXD getHitPatternVXD() const
Getter for the hit pattern in the VXD;.
This class represents an ideal helix in perigee parameterization including the covariance matrix of t...
TMatrixDSym getCartesianCovariance(const double bZ_tesla=1.5) const
Getter for the position and momentum covariance matrix.
static const double T
[tesla]
Definition: Unit.h:120
static void getField(const double *pos, double *field)
return the magnetic field at a given position.
Definition: BFieldManager.h:91
Abstract base class for different kinds of events.
Definition: ClusterUtils.h:23