Belle II Software development
BaseFitObject.h
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * Forked from https://github.com/iLCSoft/MarlinKinfit *
6 * *
7 * Further information about the fit engine and the user interface *
8 * provided in MarlinKinfit can be found at *
9 * https://www.desy.de/~blist/kinfit/doc/html/ *
10 * and in the LCNotes LC-TOOL-2009-001 and LC-TOOL-2009-004 available *
11 * from http://www-flc.desy.de/lcnotes/ *
12 * *
13 * See git log for contributors and copyright holders. *
14 * This file is licensed under LGPL-3.0, see LICENSE.md. *
15 **************************************************************************/
16
17#ifndef __BASEFITOBJECT_H
18#define __BASEFITOBJECT_H
19
20#include <iostream>
21
22#include "analysis/OrcaKinFit/BaseDefs.h"
23
24// Class BaseFitObject
26
83namespace Belle2 {
89 namespace OrcaKinFit {
90
92 public:
95
98 );
101 ) ;
102
104 virtual ~BaseFitObject();
105
107 virtual BaseFitObject* copy() const = 0;
108
110 virtual BaseFitObject& assign(const BaseFitObject& source
111 );
112
114 virtual bool setParam(int ilocal,
115 double par_,
116 bool measured_,
117 bool fixed_ = false
118 );
119
121 virtual bool setParam(int ilocal,
122 double par_
123 );
124
126 virtual bool updateParams(double p[],
127 int idim
128 );
129
131 virtual bool setMParam(int ilocal,
132 double mpar_
133 );
134
136 virtual bool setError(int ilocal,
137 double err_
138 );
139
141 virtual bool setCov(int ilocal,
142 int jlocal,
143 double cov_
144 );
145
148 virtual bool setGlobalParNum(int ilocal,
149 int iglobal
150 );
151
153 virtual bool fixParam(int ilocal,
154 bool fix = true
155 );
157 virtual bool releaseParam(int ilocal
158 )
159 { return fixParam(ilocal, false); }
160
162 virtual bool isParamFixed(int ilocal
163 ) const;
164
166 virtual double getParam(int ilocal
167 ) const;
169 virtual const char* getParamName(int ilocal
170 ) const = 0;
172 // virtual const char *getName () const { return name ? name : "???";}
173 virtual const char* getName() const; // { return name ? name : "???";}
175 void setName(const char* name_);
177 virtual double getMParam(int ilocal
178 ) const;
180 virtual double getError(int ilocal
181 ) const;
183 virtual double getCov(int ilocal,
184 int jlocal
185 ) const;
187 virtual double getRho(int ilocal,
188 int jlocal
189 ) const;
191 virtual bool isParamMeasured(int ilocal
192 ) const;
194 virtual int getGlobalParNum(int ilocal
195 ) const;
197 virtual int getNPar() const = 0;
199 virtual int getNMeasured() const;
201 virtual int getNUnmeasured() const;
203 virtual int getNFree() const;
205 virtual int getNFixed() const;
206
208 virtual double getChi2() const;
209
211 virtual double getDChi2DParam(int ilocal
212 ) const ;
213
215 virtual double getD2Chi2DParam2(int ilocal,
216 int jlocal
217 ) const;
218
220 virtual std::ostream& printParams(std::ostream& os
221 ) const;
222
224 virtual std::ostream& printRhoValues(std::ostream& os
225 ) const;
226
228 virtual std::ostream& print1stDerivatives(std::ostream& os
229 ) const;
230
232 virtual std::ostream& print2ndDerivatives(std::ostream& os
233 ) const;
234
236 virtual std::ostream& print(std::ostream& os
237 ) const = 0;
238
239
241 void invalidateCache() const {cachevalid = false;};
242 virtual void updateCache() const = 0;
243
244 // these are the methods that fill the fitter's matrices/vectors
245
248 virtual void addToGlobCov(double* glcov,
249 int idim
250 ) const;
251
253 virtual int addToGlobalChi2DerVector(double* y,
254 int idim
255 ) const;
256
258 virtual void addToGlobalChi2DerMatrix(double* M,
259 int idim
260 ) const;
261
263 virtual void addToGlobalChi2DerVector(double* y,
264 int idim,
265 double lambda,
266 double der[],
267 int metaSet
268 ) const;
269
270 // seems not to be used DJeans
271 // virtual void addToDerivatives (double der[], int idim,
272 // double factor[], int metaSet) const;
273
274 virtual void addTo1stDerivatives(double M[],
275 int idim,
276 double der[],
277 int kglobal,
278 int metaSet
279 ) const;
280
281 virtual void addTo2ndDerivatives(double der2[], int idim, double factor[], int metaSet) const;
282 virtual void addTo2ndDerivatives(double M[], int idim, double lambda, double der[], int metaSet) const;
283
284 // DANIEL added
285 // derivatives of intermediate variables wrt object's local parameters
286 // these must be implemented by the derived classes for each type of object
287 virtual double getFirstDerivative_Meta_Local(int iMeta,
288 int ilocal,
289 int metaSet
290 ) const = 0;
291
292 virtual double getSecondDerivative_Meta_Local(int iMeta,
293 int ilocal,
294 int jlocal,
295 int metaSet
296 ) const = 0;
297
298 virtual void initCov();
299
300 virtual double getError2(double der[], int metaset) const;
301
302 virtual void getDerivatives(double der[], int idim) const = 0;
303
304
305 protected:
306 char* name;
307 const static double eps2;
308
309 // DANIEL moved all of this stuff to BaseFitObject
310 // it costs some extra memory, since everything has dimension of largest # of parameters
311 // but avoid a lot of almost-duplication in the derived classes
312
314 virtual bool calculateCovInv() const;
315
317 double par[BaseDefs::MAXPAR];
319 double mpar[BaseDefs::MAXPAR];
321 bool measured[BaseDefs::MAXPAR];
323 bool fixed[BaseDefs::MAXPAR];
325 int globalParNum [BaseDefs::MAXPAR];
327 double cov [BaseDefs::MAXPAR][BaseDefs::MAXPAR];
329 mutable double covinv [BaseDefs::MAXPAR][BaseDefs::MAXPAR];
331 mutable bool covinvvalid;
333 mutable bool cachevalid;
334 // end DANIEL adds
335
336 };
337
341 inline std::ostream& operator<< (std::ostream& os,
342 const BaseFitObject& bfo
343 )
344 {
345 return bfo.print(os);
346 }
347
348 }// end OrcaKinFit namespace
350} // end Belle2 namespace
351
352
353#endif // __BASEFITOBJECT_H
354
bool covinvvalid
flag for valid inverse covariance matrix
virtual ~BaseFitObject()
Virtual destructor.
virtual bool updateParams(double p[], int idim)
Read values from global vector, readjust vector; return: significant change.
BaseFitObject & operator=(const BaseFitObject &rhs)
Assignment.
virtual int getNPar() const =0
Get total number of parameters of this FitObject.
virtual void addToGlobalChi2DerMatrix(double *M, int idim) const
Add 2nd derivatives of chi squared to global derivative matrix.
virtual const char * getParamName(int ilocal) const =0
Get name of parameter ilocal.
BaseFitObject()
Default constructor.
bool fixed[BaseDefs::MAXPAR]
fixed flag
virtual BaseFitObject & assign(const BaseFitObject &source)
Assign from anther object, if of same type.
virtual double getError(int ilocal) const
Get error of parameter ilocal.
virtual double getRho(int ilocal, int jlocal) const
Get correlation coefficient between parameters ilocal and jlocal.
virtual std::ostream & printRhoValues(std::ostream &os) const
print the correlation coefficients
virtual int getNFixed() const
Get number of fixed parameters of this FitObject.
virtual bool setParam(int ilocal, double par_, bool measured_, bool fixed_=false)
Set value and measured flag of parameter i; return: significant change.
void setName(const char *name_)
Set object's name.
double covinv[BaseDefs::MAXPAR][BaseDefs::MAXPAR]
inverse pf local covariance matrix
virtual double getCov(int ilocal, int jlocal) const
Get covariance between parameters ilocal and jlocal.
virtual double getChi2() const
Get chi squared from measured and fitted parameters.
virtual double getD2Chi2DParam2(int ilocal, int jlocal) const
Get second derivative of chi squared w.r.t. parameters ilocal1 and ilocal2.
bool cachevalid
flag for valid cache
virtual bool fixParam(int ilocal, bool fix=true)
Fix a parameter (fix=true), or release it (fix=false)
int globalParNum[BaseDefs::MAXPAR]
global parameter number for each parameter
virtual bool setCov(int ilocal, int jlocal, double cov_)
Set covariance of parameters ilocal and jlocal; return: success.
virtual std::ostream & print(std::ostream &os) const =0
print object to ostream
virtual int getGlobalParNum(int ilocal) const
Get global parameter number of parameter ilocal.
virtual bool setMParam(int ilocal, double mpar_)
Set measured value of parameter ilocal; return: success.
virtual const char * getName() const
Get object's name.
virtual int getNUnmeasured() const
Get number of unmeasured parameters of this FitObject.
virtual std::ostream & print1stDerivatives(std::ostream &os) const
print the 1st derivatives wrt metaSet 0 (E, px, py, pz)
virtual bool isParamFixed(int ilocal) const
Returns whether parameter is fixed.
double par[BaseDefs::MAXPAR]
fit parameters
virtual std::ostream & print2ndDerivatives(std::ostream &os) const
print the 2nd derivatives wrt metaSet 0 (E, px, py, pz)
void invalidateCache() const
invalidate any cached quantities
double mpar[BaseDefs::MAXPAR]
measured parameters
virtual double getParam(int ilocal) const
Get current value of parameter ilocal.
double cov[BaseDefs::MAXPAR][BaseDefs::MAXPAR]
local covariance matrix
bool measured[BaseDefs::MAXPAR]
measured flag
virtual int addToGlobalChi2DerVector(double *y, int idim) const
Add derivatives of chi squared to global derivative vector.
virtual void addToGlobCov(double *glcov, int idim) const
Add covariance matrix elements to global covariance matrix of size idim x idim.
virtual std::ostream & printParams(std::ostream &os) const
print the parameters and errors
virtual bool calculateCovInv() const
Calculate the inverse of the covariance matrix.
virtual bool isParamMeasured(int ilocal) const
Get measured flag for parameter ilocal.
virtual bool setGlobalParNum(int ilocal, int iglobal)
Set number of parameter ilocal in global list return true signals OK.
virtual bool releaseParam(int ilocal)
Release a parameter.
virtual double getMParam(int ilocal) const
Get measured value of parameter ilocal.
virtual double getDChi2DParam(int ilocal) const
Get derivative of chi squared w.r.t. parameter ilocal.
virtual int getNFree() const
Get number of free parameters of this FitObject.
virtual int getNMeasured() const
Get number of measured parameters of this FitObject.
virtual bool setError(int ilocal, double err_)
Set error of parameter ilocal; return: success.
virtual BaseFitObject * copy() const =0
Return a new copy of itself.
Abstract base class for different kinds of events.