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

Last change on this file since 396 was 307, checked in by ansari, 26 years ago

FIO_... + grosses modifs cmv 19/5/99
File size: 7.8 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 "matrix.h"
7#include "cvector.h"
8#include "generaldata.h"
9
10namespace PlanckDPC {
11
12//================================================================
13// GeneralFunction
14//================================================================
15
16class GeneralFunction {
17public:
18 GeneralFunction(unsigned int nVar, unsigned int nPar);
19 virtual ~GeneralFunction();
20
21 virtual double Value(double const xp[], double const* parm)=0;
22 virtual double Val_Der(double const xp[], double const* parm
23 , double* DgDpar);
24
25 void SetDeltaParm(int numPar, double delta=0.);
26 void SetDeltaParm(double const* dparam);
27
28 inline int NVar() const {return mNVar;}
29 inline int NPar() const {return mNPar;}
30
31protected:
32 const int mNVar; // nombre de variables f(x,y,z,...)
33 const int mNPar; // nombre de parametres
34
35 double *deltaParm;
36 double *tmpParm;
37};
38
39//================================================================
40// GeneralFunc
41//================================================================
42
43class GeneralFunc : public GeneralFunction {
44public:
45 GeneralFunc(unsigned int nvar, unsigned int npar
46 ,double (*fun)(double const*,double const*)
47 ,double (*funder)(double const*, double const*, double*)=NULL);
48 virtual ~GeneralFunc();
49
50 virtual double Value(double const xp[], double const* Par);
51 virtual double Val_Der(double const xp[],double const* parm
52 , double* DgDpar);
53
54protected:
55 double (*tmpFun) (double const*, double const*);
56 double (*tmpFunDer)(double const*, double const*, double*);
57};
58
59//================================================================
60// GeneralXi2
61//================================================================
62
63class GeneralFitData;
64
65class GeneralXi2 {
66public:
67 GeneralXi2(unsigned int nPar);
68 virtual ~GeneralXi2();
69
70 virtual double Value(GeneralFitData& data, double* parm, int& ndataused)=0;
71 virtual double Derivee(GeneralFitData& data, int i, double* parm);
72 virtual double Derivee2(GeneralFitData& data, int i,int j, double* parm);
73
74 void SetDeltaParm(int numPar, double delta=0.);
75 void SetDeltaParm(double const* dparam);
76
77 inline int NPar() const {return mNPar;}
78
79protected:
80 const int mNPar; // nombre de parametres
81
82 double *deltaParm;
83};
84
85//================================================================
86// GENERALFIT
87//================================================================
88
89class GeneralFit {
90public:
91 GeneralFit(GeneralFunction* f);
92 GeneralFit(GeneralXi2* f);
93 ~GeneralFit();
94
95 void WriteStep(char *filename = NULL);
96 void SetDebug(int level = 0);
97 void SetMaxStep(int n = 100);
98 void SetLambda_Fac(double fac = 10.);
99 void SetStopChi2(double s = 0.01);
100 void SetEps(double ep = 1.e-8);
101 void SetEps(int n,double ep = 1.e-8);
102 void SetStopMx(int nstopmx = 3, double stopchi2 = -1.);
103 void SetStopLent(int nstoplent = 3);
104 void SetFunction(GeneralFunction*);
105 void SetFuncXi2(GeneralXi2*);
106 void SetData(GeneralFitData*);
107 void SetParam(int n,double value, double step
108 ,double min=1., double max=-1.);
109 void SetParam(int n,string const&
110 ,double value, double step
111 ,double min=1., double max=-1.);
112 void SetParam(int n,double value);
113 void SetStep(int n,double value);
114 void SetMinStepDeriv(int i,double val = 0.);
115 void SetMinStepDeriv(double val = 0.);
116 void SetBound(int n,double min,double max);
117 void SetBound(int n);
118 void SetUnBound(int n);
119 void SetUnBound();
120 void SetFix(int n,double v);
121 void SetFix(int n);
122 void SetFree(int n);
123 void SetFree();
124
125 double GetParm(int);
126 Vector GetParm();
127 double GetParmErr(int);
128 double GetCoVar(int,int);
129 double GetStep(int n);
130 double GetMax(int n);
131 double GetMin(int n);
132 inline double GetChi2() const {return Chi2;};
133 inline double GetChi2Red() const {
134 if(mNddl<=0) return (double) mNddl;
135 return Chi2/(double) mNddl;
136 };
137 inline double GetEps(int i) const {return Eps(i);};
138 inline int GetNddl() const {return mNddl;};
139 inline int GetNStep() const {return nStep;};
140 inline int GetNVar() const {return mNVar;};
141 inline int GetNPar() const {return mNPar;};
142 inline int GetNFree() const {return mNParFree;};
143 inline int GetNBound() const {return mNParBound;};
144 inline int GetNStop() const {return nStop;};
145 inline int GetNStopLent() const {return nStopLent;};
146 inline GeneralFunction* GetFunction() const {return mFunction;};
147 inline GeneralFitData* GetGData() const {return mData;};
148
149 void PrintStatus();
150 void PrintFit();
151 void PrintParm(int n);
152 void PrintParm();
153
154 int Fit();
155 double ReCalChi2(int& nddl, double* par = NULL);
156 GeneralFitData DataResidus(bool clean=true);
157 GeneralFitData DataFunction(bool clean=true);
158 void PrintFitErr(int rc);
159
160protected:
161 int mNtry; // numero d'appel de la routine de fit.
162 int mNVar; // nombre de variables f(x,y,z,...)
163 int mNPar; // nombre de parametres
164 int mNParFree; // nombre de parametres libres
165 int mNParBound; // nombre de parametres bornes
166 GeneralFunction* mFunction;
167 GeneralXi2* mFuncXi2;
168 GeneralFitData* mData;
169
170 Vector Param;
171 Vector errParam;
172 Vector stepParam;
173 Vector minParam;
174 Vector maxParam;
175 Vector minStepDeriv;
176 Vector Eps;
177 unsigned short int* fixParam;
178 unsigned short int* boundParam;
179 string* nameParam;
180
181 double Lambda_Fac;
182 double stopChi2;
183 int maxStep;
184 int nStopMx;
185 double stopChi2SMx;
186 int nStopLent;
187 int debugLevel;
188 FILE *FileStep;
189
190 Matrix ATGA;
191 Vector BETA;
192 Matrix ATGA_Try;
193 Vector BETA_Try;
194 Vector C;
195 Vector D;
196
197 double Chi2;
198 int mNddl;
199 int nStep;
200 int nStop, nStopL;
201 double Lambda;
202
203 // Fonctions privees
204 void write_in_step(double ci2,Vector& par);
205 void General_Init(void);
206 void TryFunc(Vector& par,Vector& par_tr);
207 void TryXi2(Vector& par,Vector& par_tr);
208 void CheckSanity();
209 void Set_Bound_C_D(int i);
210 void Set_Bound_C_D();
211 double p_vers_tr(int i,double p);
212 Vector p_vers_tr(Vector const& p);
213 void p_vers_tr(Vector const& p,Vector& tr);
214 double tr_vers_p(int i,double tr);
215 Vector tr_vers_p(Vector const& tr);
216 void tr_vers_p(Vector const& tr,Vector& p);
217 double c_dp_vers_dtr(int i,double tr);
218 Vector dp_vers_dtr(Vector const& dp,Vector const& tr);
219 void dp_vers_dtr(Vector const& dp,Vector const& tr,Vector& dtr);
220 double c_dtr_vers_dp(int i,double tr);
221 Vector dtr_vers_dp(Vector const& dtr,Vector const& tr);
222 int put_in_limits_for_deriv(Vector const& p,Vector& dp,double dist=0.66);
223 inline void dtr_vers_dp(Vector const& dtr,Vector const& tr,Vector& dp)
224 { for(int i=0;i<mNPar;i++)
225 { if( fixParam[i] ) continue;
226 if( ! boundParam[i] ) dp(i) = dtr(i);
227 else dp(i) = D(i)/(1.+tr(i)*tr(i)) * dtr(i); }
228 };
229};
230
231} // Fin du namespace
232
233#endif
Note: See TracBrowser for help on using the repository browser.