source: Sophya/trunk/Poubelle/DPC:FitsIOServer/NTools/generalfit.h@ 3038

Last change on this file since 3038 was 658, checked in by ansari, 26 years ago

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