source: Sophya/trunk/SophyaLib/NTools/generalfit.h@ 914

Last change on this file since 914 was 914, checked in by ansari, 25 years ago

documentation cmv 13/4/00
File size: 9.4 KB
Line 
1// This may look like C code, but it is really -*- C++ -*-
2#ifndef GENERALFIT_SEEN
3#define GENERALFIT_SEEN
4
5#include "pexceptions.h"
6#include "tvector.h"
7#include "generaldata.h"
8
9namespace SOPHYA {
10
11//================================================================
12// GeneralFunction
13//================================================================
14
15
16/*!
17 Classe de fonctions parametrees a plusieurs variables:
18 \f$ F[x1,x2,x3,...:a1,a2,a3,...] \f$
19*/
20class GeneralFunction {
21public:
22 GeneralFunction(unsigned int nVar, unsigned int nPar);
23 virtual ~GeneralFunction();
24
25 //! Valeur de la fonction a definir par l'utilisateur (fct virtuelle pure)
26 virtual double Value(double const xp[], double const* parm)=0;
27 virtual double Val_Der(double const xp[], double const* parm
28 , double* DgDpar);
29
30 void SetDeltaParm(int numPar, double delta=0.);
31 void SetDeltaParm(double const* dparam);
32
33 //! Retourne le nombre de variables Xi
34 inline int NVar() const {return mNVar;}
35 //! Retourne le nombre de parametres Ai
36 inline int NPar() const {return mNPar;}
37
38protected:
39 const int mNVar; //!< nombre de variables f(x,y,z,...)
40 const int mNPar; //!< nombre de parametres
41
42 double *deltaParm;
43 double *tmpParm;
44};
45
46//================================================================
47// GeneralFunc
48//================================================================
49
50/*!
51 Classe de fonctions parametrees a plusieurs variables
52 derivant de ``GeneralFunction''. Permet de definir
53 une fonction a fiter sans passer par une classe derivee
54 en utilisant l'ecriture courante du C. La fonction
55 retournant les derivees par rapport aux parametres du fit
56 peut etre egalement fournie (optionnel).
57*/
58class GeneralFunc : public GeneralFunction {
59public:
60 GeneralFunc(unsigned int nvar, unsigned int npar
61 ,double (*fun)(double const*,double const*)
62 ,double (*funder)(double const*, double const*, double*)=NULL);
63 virtual ~GeneralFunc();
64
65 virtual double Value(double const xp[], double const* Par);
66 virtual double Val_Der(double const xp[],double const* parm
67 , double* DgDpar);
68
69protected:
70 double (*tmpFun) (double const*, double const*);
71 double (*tmpFunDer)(double const*, double const*, double*);
72};
73
74//================================================================
75// GeneralXi2
76//================================================================
77
78class GeneralFitData;
79
80/*!
81 Classe de Xi2 a plusieurs parametres :
82 \f$ Xi2[a1,a2,a3,...] \f$
83*/
84class GeneralXi2 {
85public:
86 GeneralXi2(unsigned int nPar);
87 virtual ~GeneralXi2();
88
89 /*!
90 Valeur du Xi2 a definir par l'utilisateur (fct virtuelle pure)
91 a partir des donnees de `data'. l'utilisateur doit egalement
92 retourner le nombre de points de mesure utilises dans le calcul
93 du Xi2 (`ndataused').
94 */
95 virtual double Value(GeneralFitData& data, double* parm, int& ndataused)=0;
96 virtual double Derivee(GeneralFitData& data, int i, double* parm);
97 virtual double Derivee2(GeneralFitData& data, int i,int j, double* parm);
98
99 void SetDeltaParm(int numPar, double delta=0.);
100 void SetDeltaParm(double const* dparam);
101
102 //! Retourne le nombre de parametres Ai.
103 inline int NPar() const {return mNPar;}
104
105protected:
106 const int mNPar; //!< nombre de parametres
107
108 double *deltaParm;
109};
110
111//================================================================
112// GENERALFIT
113//================================================================
114
115//! Classe de fit d'une GeneralFunction sur une GeneralFitData
116class GeneralFit {
117public:
118 GeneralFit(GeneralFunction* f);
119 GeneralFit(GeneralXi2* f);
120 ~GeneralFit();
121
122 void WriteStep(char *filename = NULL);
123 void SetDebug(int level = 0);
124 void SetMaxStep(int n = 100);
125 void SetLambda_Fac(double fac = 10.);
126 void SetStopChi2(double s = 0.01);
127 void SetEps(double ep = 1.e-8);
128 void SetEps(int n,double ep = 1.e-8);
129 void SetStopMx(int nstopmx = 3, double stopchi2 = -1.);
130 void SetStopLent(int nstoplent = 3);
131 void SetFunction(GeneralFunction*);
132 void SetFuncXi2(GeneralXi2*);
133 void SetData(GeneralFitData*);
134 void SetParam(int n,double value, double step
135 ,double min=1., double max=-1.);
136 void SetParam(int n,string const&
137 ,double value, double step
138 ,double min=1., double max=-1.);
139 void SetParam(int n,double value);
140 void SetStep(int n,double value);
141 void SetMinStepDeriv(int i,double val = 0.);
142 void SetMinStepDeriv(double val = 0.);
143 void SetBound(int n,double min,double max);
144 void SetBound(int n);
145 void SetUnBound(int n);
146 void SetUnBound();
147 void SetFix(int n,double v);
148 void SetFix(int n);
149 void SetFree(int n);
150 void SetFree();
151
152 double GetParm(int);
153 Vector GetParm();
154 double GetParmErr(int);
155 double GetCoVar(int,int);
156 double GetStep(int n);
157 double GetMax(int n);
158 double GetMin(int n);
159 //! Retourne le Chi2
160 inline double GetChi2() const {return Chi2;};
161 //! Retourne le Chi2 reduit
162 inline double GetChi2Red() const {
163 if(mNddl<=0) return (double) mNddl;
164 return Chi2/(double) mNddl;
165 };
166 //! Retourne la precision de convergence pour le parametre i.
167 inline double GetEps(int i) const {return Eps(i);};
168 //! Retourne le nombre de degres de liberte
169 inline int GetNddl() const {return mNddl;};
170 //! Retourne le nombre d'iterations
171 inline int GetNStep() const {return nStep;};
172 //! Retourne le nombre de variables
173 inline int GetNVar() const {return mNVar;};
174 //! Retourne le nombre de parametres
175 inline int GetNPar() const {return mNPar;};
176 //! Retourne le nombre de parametres libres
177 inline int GetNFree() const {return mNParFree;};
178 //! Retourne le nombre de parametres bornes
179 inline int GetNBound() const {return mNParBound;};
180 //! Retourne le nstop de convergence
181 inline int GetNStop() const {return nStop;};
182 //! Retourne le nstop de convergence lente.
183 inline int GetNStopLent() const {return nStopLent;};
184 //! Retourne le pointeur sur la GeneralFunction utilisee.
185 inline GeneralFunction* GetFunction() const {return mFunction;};
186 //! Retourne le pointeur sur la GeneralFitData utilisee.
187 inline GeneralFitData* GetGData() const {return mData;};
188
189 void PrintStatus();
190 void PrintFit();
191 void PrintParm(int n);
192 void PrintParm();
193
194 int Fit();
195 double ReCalChi2(int& nddl, double* par = NULL);
196 GeneralFitData DataResidus(bool clean=true);
197 GeneralFitData DataFunction(bool clean=true);
198 void PrintFitErr(int rc);
199
200protected:
201 int mNtry; //!< numero d'appel de la routine de fit.
202 int mNVar; //!< nombre de variables f(x,y,z,...)
203 int mNPar; //!< nombre de parametres
204 int mNParFree; //!< nombre de parametres libres
205 int mNParBound; //!< nombre de parametres bornes
206 GeneralFunction* mFunction;
207 GeneralXi2* mFuncXi2;
208 GeneralFitData* mData;
209
210 Vector Param;
211 Vector errParam;
212 Vector stepParam;
213 Vector minParam;
214 Vector maxParam;
215 Vector minStepDeriv;
216 Vector Eps;
217 unsigned short int* fixParam;
218 unsigned short int* boundParam;
219 string* nameParam;
220
221 double Lambda_Fac;
222 double stopChi2;
223 int maxStep;
224 int nStopMx;
225 double stopChi2SMx;
226 int nStopLent;
227 int debugLevel;
228 FILE *FileStep;
229
230 Matrix ATGA;
231 Vector BETA;
232 Matrix ATGA_Try;
233 Vector BETA_Try;
234 Vector C;
235 Vector D;
236
237 double Chi2;
238 int mNddl;
239 int nStep;
240 int nStop, nStopL;
241 double Lambda;
242
243 // Fonctions privees
244 void write_in_step(double ci2,Vector& par);
245 void General_Init(void);
246 void TryFunc(Vector& par,Vector& par_tr);
247 void TryXi2(Vector& par,Vector& par_tr);
248 void CheckSanity();
249 void Set_Bound_C_D(int i);
250 void Set_Bound_C_D();
251 double p_vers_tr(int i,double p);
252 Vector p_vers_tr(Vector const& p);
253 void p_vers_tr(Vector const& p,Vector& tr);
254 double tr_vers_p(int i,double tr);
255 Vector tr_vers_p(Vector const& tr);
256 void tr_vers_p(Vector const& tr,Vector& p);
257 double c_dp_vers_dtr(int i,double tr);
258 Vector dp_vers_dtr(Vector const& dp,Vector const& tr);
259 void dp_vers_dtr(Vector const& dp,Vector const& tr,Vector& dtr);
260 double c_dtr_vers_dp(int i,double tr);
261 Vector dtr_vers_dp(Vector const& dtr,Vector const& tr);
262 int put_in_limits_for_deriv(Vector const& p,Vector& dp,double dist=0.66);
263 inline void dtr_vers_dp(Vector const& dtr,Vector const& tr,Vector& dp)
264 { for(int i=0;i<mNPar;i++)
265 { if( fixParam[i] ) continue;
266 if( ! boundParam[i] ) dp(i) = dtr(i);
267 else dp(i) = D(i)/(1.+tr(i)*tr(i)) * dtr(i); }
268 };
269};
270
271} // Fin du namespace
272
273#endif
Note: See TracBrowser for help on using the repository browser.