Belle II Software development
KalmanCalculator.h
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * External Contributor: Wouter Hulsbergen *
5 * *
6 * See git log for contributors and copyright holders. *
7 * This file is licensed under LGPL-3.0, see LICENSE.md. *
8 **************************************************************************/
9#pragma once
10
11#include <analysis/VertexFitting/TreeFitter/EigenStackConfig.h>
12#include <analysis/VertexFitting/TreeFitter/ErrCode.h>
13
14// inverse() needs this, in the other classes we get away with just Eigen/Core
15#include <Eigen/Dense>
16
17
18namespace TreeFitter {
19 class FitParams;
20
23 public:
24
27 int sizeRes,
28 int sizeState
29 );
30
32 void updateState(FitParams& fitparams);
33
35 void updateState(FitParams& fitparams, FitParams& oldState);
36
38 void updateCovariance(FitParams& fitparams);
39
41 double getChiSquare() const { return m_chisq;}
42
45 const Eigen::Matrix < double, -1, 1, 0, 7, 1 > & residuals,
46 const Eigen::Matrix < double, -1, -1, 0, 7, MAX_MATRIX_SIZE > & G,
47 const FitParams& fitparams,
48 const Eigen::Matrix < double, -1, -1, 0, 7, 7 > * V = 0,
49 double weight = 1);
50
52 double getConstraintDim() const { return m_constrDim; }
53
54 private:
57
59 double m_chisq;
60
71
73 Eigen::Matrix < double, -1, 1, 0, 7, 1 > m_res;
74
76 Eigen::Matrix < double, -1, -1, 0, 7, MAX_MATRIX_SIZE > m_G;
77
79 Eigen::Matrix < double, -1, -1, 0, 7, 7 > m_R;
80
82 Eigen::Matrix < double, -1, -1, 0, 7, 7 > m_Rinverse;
83
85 Eigen::Matrix < double, -1, -1, 0, MAX_MATRIX_SIZE, 7 > m_K;
86
88 Eigen::Matrix < double, -1, -1, 0, MAX_MATRIX_SIZE, 7 > m_CGt;
89
90 };
91}
abstract errorocode be aware that the default is success
Definition ErrCode.h:14
Class to store and manage fitparams (statevector)
Definition FitParams.h:20
void updateState(FitParams &fitparams)
update statevector
ErrCode calculateGainMatrix(const Eigen::Matrix< double, -1, 1, 0, 7, 1 > &residuals, const Eigen::Matrix< double, -1, -1, 0, 7, MAX_MATRIX_SIZE > &G, const FitParams &fitparams, const Eigen::Matrix< double, -1, -1, 0, 7, 7 > *V=0, double weight=1)
init the kalman machienery
Eigen::Matrix< double, -1, -1, 0, 7, 7 > m_Rinverse
R inverse.
double getConstraintDim() const
get dimension of the constraint
KalmanCalculator(int sizeRes, int sizeState)
constructor
Eigen::Matrix< double, -1, -1, 0, 7, 7 > m_R
R residual covariance.
Eigen::Matrix< double, -1, -1, 0, 7, MAX_MATRIX_SIZE > m_G
G former H, transforms covraince of {residuals}<->{x,p,E}.
Eigen::Matrix< double, -1, 1, 0, 7, 1 > m_res
we know the max sizes of the matrices we assume the tree is smaller than MAX_MATRIX_SIZE parameters a...
void updateCovariance(FitParams &fitparams)
update the statevectors covariance
Eigen::Matrix< double, -1, -1, 0, MAX_MATRIX_SIZE, 7 > m_K
K kalman gain matrix.
double getChiSquare() const
get chi2 of this iteration
int m_constrDim
dimension of the constraint
Eigen::Matrix< double, -1, -1, 0, MAX_MATRIX_SIZE, 7 > m_CGt
C times G^t.