MarlinKinfit
0.4.0
|
A kinematic fitter using the Newton-Raphson method to solve the equations. More...
#include <NewFitterGSL.h>
Public Types | |
enum | { NPARMAX =50, NCONMAX =10, NUNMMAX =10 } |
enum | { NITMAX = 100 } |
Public Member Functions | |
NewFitterGSL () | |
Constructor. | |
virtual | ~NewFitterGSL () |
Virtual destructor. | |
virtual double | fit () |
The fit method, returns the fit probability. | |
virtual int | getError () const |
Get the error code of the last fit: 0=OK, 1=failed. | |
virtual double | getProbability () const |
Get the fit probability of the last fit. | |
virtual double | getChi2 () const |
Get the chi**2 of the last fit. | |
virtual int | getDoF () const |
Get the number of degrees of freedom of the last fit. | |
virtual int | getIterations () const |
Get the number of iterations of the last fit. | |
virtual int | getNcon () const |
Get the number of hard constraints of the last fit. | |
virtual int | getNsoft () const |
Get the number of soft constraints of the last fit. | |
virtual int | getNpar () const |
Get the number of all parameters of the last fit. | |
virtual int | getNunm () const |
Get the number of unmeasured parameters of the last fit. | |
virtual bool | initialize () |
Initialize the fitter. | |
virtual void | setDebug (int debuglevel) |
Set the Debug Level. | |
virtual void | determineLambdas (gsl_vector *vecxnew, const gsl_matrix *MatM, const gsl_vector *vecx, gsl_matrix *MatW, gsl_vector *vecw, double eps=0) |
Determine best lambda values. More... | |
virtual void | calc2ndOrderCorr (gsl_vector *vecdxhat, const gsl_vector *vecxnew, const gsl_matrix *MatM, gsl_matrix *MatW, gsl_vector *vecw, double eps=0) |
Calculate 2nd order correction step. More... | |
virtual gsl_matrix_view | calcZ (int &rankA, gsl_matrix *MatW1, gsl_matrix *MatW2, gsl_vector *vecw1, gsl_vector *vecw2, gsl_permutation *permW, double eps=0) |
Calculate null space of constraints; return value is a matrix view of MatW1, with column vectors spanning null(A) More... | |
virtual gsl_matrix_view | calcReducedHessian (int &rankH, gsl_matrix *MatW1, const gsl_vector *vecx, gsl_matrix *MatW2, gsl_matrix *MatW3, gsl_vector *vecw1, gsl_vector *vecw2, gsl_permutation *permW, double eps=0) |
Calculate reduced Hessian; return value is a matrix view of MatW1 giving the reduced Hessian. More... | |
virtual gsl_vector_view | calcReducedHessianEigenvalues (int &rankH, gsl_matrix *MatW1, const gsl_vector *vecx, gsl_matrix *MatW2, gsl_matrix *MatW3, gsl_vector *vecw1, gsl_vector *vecw2, gsl_permutation *permW, gsl_eigen_symm_workspace *eigenws, double eps=0) |
virtual double | calcChi2 () |
Calculate the chi2. | |
void | fillx (gsl_vector *vecx) |
void | fillperr (gsl_vector *vece) |
void | assembleM (gsl_matrix *MatM, const gsl_vector *vecx, bool errorpropagation=false) |
void | assembleG (gsl_matrix *MatM, const gsl_vector *vecx) |
void | scaleM (gsl_matrix *MatMscal, const gsl_matrix *MatM, const gsl_vector *vece) |
void | assembley (gsl_vector *vecy, const gsl_vector *vecx) |
void | scaley (gsl_vector *vecyscal, const gsl_vector *vecy, const gsl_vector *vece) |
void | assembleChi2Der (gsl_vector *vecy) |
void | addConstraints (gsl_vector *vecy) |
void | assembleConstDer (gsl_matrix *MatM) |
int | calcNewtonDx (gsl_vector *vecdx, gsl_vector *vecdxscal, gsl_vector *vecx, const gsl_vector *vece, gsl_matrix *MatM, gsl_matrix *MatMscal, gsl_vector *vecy, gsl_vector *vecyscal, gsl_matrix *MatW, gsl_matrix *MatW2, gsl_permutation *permW, gsl_vector *vecw) |
int | calcLimitedDx (double &alpha, double &mu, gsl_vector *vecxnew, int imode, gsl_vector *vecx, gsl_vector *vecdxhat, const gsl_vector *vecdx, const gsl_vector *vecdxscal, const gsl_vector *vece, const gsl_matrix *MatM, const gsl_matrix *MatMscal, gsl_matrix *MatW, gsl_vector *vecw) |
int | doLineSearch (double &alpha, gsl_vector *vecxnew, int imode, double phi0, double dphi0, double phiR, double eta, double zeta, double mu, const gsl_vector *vecx, const gsl_vector *vecdx, const gsl_vector *vece, gsl_vector *vecw) |
double | calcMu (const gsl_vector *vecx, const gsl_vector *vece, const gsl_vector *vecdx, const gsl_vector *vecdxscal, const gsl_vector *xnew, const gsl_matrix *MatM, const gsl_matrix *MatMscal, gsl_vector *vecw) |
double | meritFunction (double mu, const gsl_vector *vecx, const gsl_vector *vece) |
double | meritFunctionDeriv (double mu, const gsl_vector *vecx, const gsl_vector *vece, const gsl_vector *vecdx, gsl_vector *vecw) |
bool | updateParams (gsl_vector *vecx) |
int | invertM () |
void | calcCovMatrix (gsl_matrix *MatW, gsl_permutation *permW, gsl_vector *vecx) |
double | calcpTLp (const gsl_vector *vecdx, const gsl_matrix *MatM, gsl_vector *vecw) |
int | solveSystem (gsl_vector *vecdxscal, double &detW, const gsl_vector *vecyscal, const gsl_matrix *MatMscal, gsl_matrix *MatW, gsl_matrix *MatW2, gsl_vector *vecw, double epsLU, double epsSV) |
solve system of equations Mscal*dxscal = yscal | |
int | solveSystemLU (gsl_vector *vecdxscal, double &detW, const gsl_vector *vecyscal, const gsl_matrix *MatMscal, gsl_matrix *MatW, gsl_vector *vecw, double eps) |
solve system of equations Mscal*dxscal = yscal using LU decomposition | |
int | solveSystemSVD (gsl_vector *vecdxscal, double &detW, const gsl_vector *vecyscal, const gsl_matrix *MatMscal, gsl_matrix *MatW, gsl_matrix *MatW2, gsl_vector *vecw, double eps) |
solve system of equations Mscal*dxscal = yscal using SVD decomposition | |
![]() | |
virtual void | addFitObject (BaseFitObject *fitobject_) |
virtual void | addFitObject (BaseFitObject &fitobject_) |
virtual void | addConstraint (BaseConstraint *constraint_) |
virtual void | addConstraint (BaseConstraint &constraint_) |
virtual void | addHardConstraint (BaseHardConstraint *constraint_) |
virtual void | addHardConstraint (BaseHardConstraint &constraint_) |
virtual void | addSoftConstraint (BaseSoftConstraint *constraint_) |
virtual void | addSoftConstraint (BaseSoftConstraint &constraint_) |
virtual std::vector < BaseFitObject * > * | getFitObjects () |
virtual std::vector < BaseHardConstraint * > * | getConstraints () |
virtual std::vector < BaseSoftConstraint * > * | getSoftConstraints () |
virtual void | reset () |
virtual BaseTracer * | getTracer () |
virtual const BaseTracer * | getTracer () const |
virtual void | setTracer (BaseTracer *newTracer) |
virtual void | setTracer (BaseTracer &newTracer) |
virtual const double * | getGlobalCovarianceMatrix (int &idim) const |
virtual double * | getGlobalCovarianceMatrix (int &idim) |
Static Public Member Functions | |
static void | ini_gsl_permutation (gsl_permutation *&p, unsigned int size) |
static void | ini_gsl_vector (gsl_vector *&v, int unsigned size) |
static void | ini_gsl_matrix (gsl_matrix *&m, int unsigned size1, unsigned int size2) |
static void | debug_print (const gsl_matrix *m, const char *name) |
static void | debug_print (const gsl_vector *v, const char *name) |
static void | add (gsl_vector *vecz, const gsl_vector *vecx, double a, const gsl_vector *vecy) |
static bool | isfinite (const gsl_vector *vec) |
static bool | isfinite (const gsl_matrix *mat) |
static void | MoorePenroseInverse (gsl_matrix *Ainv, gsl_matrix *A, gsl_matrix *W, gsl_vector *w, double eps=0) |
Compute the Moore-Penrose pseudo-inverse A+ of A, using SVD. More... | |
Public Attributes | |
int | npar |
total number of parameters | |
int | ncon |
total number of hard constraints | |
int | nsoft |
total number of soft constraints | |
int | nunm |
total number of unmeasured parameters | |
int | ierr |
Error status. | |
int | nit |
Number of iterations. | |
double | fitprob |
fit probability | |
double | chi2 |
final chi2 | |
unsigned int | idim |
gsl_vector * | x |
gsl_vector * | xold |
gsl_vector * | xnew |
gsl_vector * | dx |
gsl_vector * | dxscal |
gsl_vector * | y |
gsl_vector * | yscal |
gsl_vector * | perr |
gsl_vector * | v1 |
gsl_vector * | v2 |
gsl_matrix * | M |
gsl_matrix * | Mscal |
gsl_matrix * | W |
gsl_matrix * | W2 |
gsl_matrix * | W3 |
gsl_matrix * | M1 |
gsl_matrix * | M2 |
gsl_matrix * | M3 |
gsl_matrix * | M4 |
gsl_matrix * | M5 |
gsl_matrix * | CC |
gsl_matrix * | CC1 |
gsl_matrix * | CCinv |
gsl_permutation * | permW |
gsl_eigen_symm_workspace * | eigenws |
unsigned int | eigenwsdim |
double | chi2best |
double | chi2new |
double | chi2old |
double | fvalbest |
double | scale |
double | scalebest |
double | stepsize |
double | stepbest |
double | scalevals [NITMAX] |
double | fvals [NITMAX] |
int | imerit |
bool | try2ndOrderCorr |
int | debug |
![]() | |
std::map< std::string, double > | traceValues |
Additional Inherited Members | |
![]() | |
typedef std::vector < BaseFitObject * > | FitObjectContainer |
typedef std::vector < BaseHardConstraint * > | ConstraintContainer |
typedef std::vector < BaseSoftConstraint * > | SoftConstraintContainer |
typedef FitObjectContainer::iterator | FitObjectIterator |
typedef ConstraintContainer::iterator | ConstraintIterator |
typedef SoftConstraintContainer::iterator | SoftConstraintIterator |
![]() | |
BaseFitter (const BaseFitter &rhs) | |
Copy constructor disabled. | |
BaseFitter & | operator= (const BaseFitter &rhs) |
Assignment disabled. | |
![]() | |
FitObjectContainer | fitobjects |
ConstraintContainer | constraints |
SoftConstraintContainer | softconstraints |
int | covDim |
dimension of global covariance matrix | |
double * | cov |
global covariance matrix of last fit problem | |
bool | covValid |
Flag whether global covariance is valid. | |
BaseTracer * | tracer |
A kinematic fitter using the Newton-Raphson method to solve the equations.
This class implements a kinematic fitter using the Newton-Raphson method to solve the system of equations arising from the Lagrange multiplier method
Author: Benno List Last update:
by:
Changelog:
|
virtual |
Calculate 2nd order correction step.
vecdxhat | the correction step |
vecxnew | the current state vector (parameters must be set to vecxnew!) |
MatM | The matrix of the last Newton step |
MatW | work matrix |
vecw | work vector |
eps | Singular values < eps*(max(abs(s_i))) are set to 0 |
Referenced by calcLimitedDx().
int NewFitterGSL::calcLimitedDx | ( | double & | alpha, |
double & | mu, | ||
gsl_vector * | vecxnew, | ||
int | imode, | ||
gsl_vector * | vecx, | ||
gsl_vector * | vecdxhat, | ||
const gsl_vector * | vecdx, | ||
const gsl_vector * | vecdxscal, | ||
const gsl_vector * | vece, | ||
const gsl_matrix * | MatM, | ||
const gsl_matrix * | MatMscal, | ||
gsl_matrix * | MatW, | ||
gsl_vector * | vecw | ||
) |
alpha | Output value alpha |
mu | Value of mu for merit function |
vecxnew | New vector x |
imode | mode: 0=Armijo, 1=Wolfe, 2=Goldstein |
vecx | Current vector x |
vecdxhat | the 2nd order correction step, if any |
vecdx | Update vector dx |
vecdxscal | Result: Update vector dx, scaled |
vece | Error vector e |
MatM | Matrix M |
MatMscal | Matrix M, scaled |
MatW | Work matrix |
vecw | Work vector w1 |
References calc2ndOrderCorr(), calcChi2(), calcMu(), doLineSearch(), meritFunction(), meritFunctionDeriv(), and BaseTracer::substep().
Referenced by fit().
double NewFitterGSL::calcMu | ( | const gsl_vector * | vecx, |
const gsl_vector * | vece, | ||
const gsl_vector * | vecdx, | ||
const gsl_vector * | vecdxscal, | ||
const gsl_vector * | xnew, | ||
const gsl_matrix * | MatM, | ||
const gsl_matrix * | MatMscal, | ||
gsl_vector * | vecw | ||
) |
vecx | Current vector x |
vece | Current errors x |
vecdx | Current step dx |
vecdxscal | Current step dx, scaled |
xnew | New vector x |
MatM | Current matrix M |
MatMscal | Current scaled matrix M |
vecw | Work vector w |
References calcpTLp(), ncon, and npar.
Referenced by calcLimitedDx().
int NewFitterGSL::calcNewtonDx | ( | gsl_vector * | vecdx, |
gsl_vector * | vecdxscal, | ||
gsl_vector * | vecx, | ||
const gsl_vector * | vece, | ||
gsl_matrix * | MatM, | ||
gsl_matrix * | MatMscal, | ||
gsl_vector * | vecy, | ||
gsl_vector * | vecyscal, | ||
gsl_matrix * | MatW, | ||
gsl_matrix * | MatW2, | ||
gsl_permutation * | permW, | ||
gsl_vector * | vecw | ||
) |
vecdx | Result: Update vector dx |
vecdxscal | Result: Update vector dx, scaled |
vecx | Current vector x |
vece | Current ``error'' set of x |
MatM | Matrix M |
MatMscal | Matrix M, scaled |
vecy | Vector y |
vecyscal | Vector y, scaled, |
MatW | Work matrix |
MatW2 | Work matrix |
permW | Work permutation vector |
vecw | Work vector |
References calcpTLp(), determineLambdas(), ncon, npar, and solveSystem().
Referenced by fit().
double NewFitterGSL::calcpTLp | ( | const gsl_vector * | vecdx, |
const gsl_matrix * | MatM, | ||
gsl_vector * | vecw | ||
) |
vecdx | Current step dx |
MatM | Current matrix M |
vecw | Work vector w |
References npar.
Referenced by calcMu(), and calcNewtonDx().
|
virtual |
Calculate reduced Hessian; return value is a matrix view of MatW1 giving the reduced Hessian.
rankH | dimension of H |
MatW1 | work matrix, contains H at the end |
vecx | vector with current x values |
MatW2 | work matrix, contains Z at the end |
MatW3 | work matrix |
vecw1 | work vector |
vecw2 | work vector |
permW | Work permutation vector |
eps | Singular values < eps*(max(abs(s_i))) are set to 0 |
Referenced by calcReducedHessianEigenvalues().
|
virtual |
rankH | dimension of H |
MatW1 | work matrix, contains H at the end |
vecx | vector with current x values |
MatW2 | work matrix, contains Z at the end |
MatW3 | work matrix |
vecw1 | work vector |
vecw2 | work vector |
permW | Work permutation vector |
eigenws | Work space for eigenvalue calculation |
eps | Singular values < eps*(max(abs(s_i))) are set to 0 |
References calcReducedHessian().
|
virtual |
Calculate null space of constraints; return value is a matrix view of MatW1, with column vectors spanning null(A)
rankA | rank of A (= number of lin. indep. constraints) |
MatW1 | work matrix, contains Z at the end |
MatW2 | work matrix |
vecw1 | work vector |
vecw2 | work vector |
permW | Work permutation vector |
eps | Singular values < eps*(max(abs(s_i))) are set to 0 |
Referenced by calcReducedHessian().
|
virtual |
Determine best lambda values.
vecxnew | vector with new lambda values |
MatM | matrix with constraint derivatives |
vecx | vector with current x values |
MatW | work matrix |
vecw | work vector |
eps | Singular values < eps*(max(abs(s_i))) are set to 0 |
Referenced by calcNewtonDx(), and fit().
int NewFitterGSL::doLineSearch | ( | double & | alpha, |
gsl_vector * | vecxnew, | ||
int | imode, | ||
double | phi0, | ||
double | dphi0, | ||
double | phiR, | ||
double | eta, | ||
double | zeta, | ||
double | mu, | ||
const gsl_vector * | vecx, | ||
const gsl_vector * | vecdx, | ||
const gsl_vector * | vece, | ||
gsl_vector * | vecw | ||
) |
alpha | Output value alpha |
vecxnew | New vector x |
imode | mode: 0=Armijo, 1=Wolfe, 2=Goldstein |
phi0 | Merit function for alpha=0 |
dphi0 | Directional derivative of merit function for alpha=0 |
phiR | Merit function for given alpha |
eta | Constant for Armijo's rule |
zeta | Constant for Wolfe's or Goldstein's rule |
mu | Value of mu for merit function |
vecx | Current vector x |
vecdx | Update vector dx |
vece | Error vector e |
vecw | Work vector w |
References meritFunction(), meritFunctionDeriv(), nit, and BaseTracer::substep().
Referenced by calcLimitedDx().
double NewFitterGSL::meritFunction | ( | double | mu, |
const gsl_vector * | vecx, | ||
const gsl_vector * | vece | ||
) |
mu | Value of mu |
vecx | Current vector x |
vece | Current errors x |
References calcChi2(), BaseHardConstraint::getGlobalNum(), and BaseHardConstraint::getValue().
Referenced by calcLimitedDx(), and doLineSearch().
double NewFitterGSL::meritFunctionDeriv | ( | double | mu, |
const gsl_vector * | vecx, | ||
const gsl_vector * | vece, | ||
const gsl_vector * | vecdx, | ||
gsl_vector * | vecw | ||
) |
mu | Value of mu |
vecx | Current vector x |
vece | Current errors x |
vecdx | Current update vector dx |
vecw | work vector |
References BaseHardConstraint::getGlobalNum(), BaseHardConstraint::getValue(), and npar.
Referenced by calcLimitedDx(), and doLineSearch().
|
static |
Compute the Moore-Penrose pseudo-inverse A+ of A, using SVD.
Ainv | Result: m x n matrix A+ |
A | Input: n x m matrix A, n >= m (is destroyed!) |
W | Work matrix, at least m x m |
w | Work vector w, at least m |
eps | Singular values < eps*(max(abs(s_i))) are set to 0 |