source: Sophya/trunk/SophyaPI/PIext/piafitting.cc@ 386

Last change on this file since 386 was 386, checked in by ercodmgr, 26 years ago

suite fit fonction user cmv 13/8/99
File size: 46.7 KB
Line 
1//CMVBUG reste pb suivant sous Linux-egcs (OSF1-cxx c'est ok):
2//openppf cmvwppf.ppf
3//fitw h1g g2 p:4500,45,10,900
4//delobjs h1g
5//clic sur FIT -> OK
6//clic sur FIT 2sd fois -> core dump!!!
7//
8//
9#include <stdio.h>
10#include <stdlib.h>
11#include <ctype.h>
12
13#include <iostream.h>
14
15#include <typeinfo>
16
17#include "piafitting.h"
18#include "strutil.h"
19
20#include "generalfit.h"
21
22#include "fct1dfit.h"
23#include "fct2dfit.h"
24
25#include "matrix.h"
26#include "cvector.h"
27#include "ntuple.h"
28#include "cimage.h"
29
30#include "histos.h"
31#include "histos2.h"
32#include "ntuple.h"
33#include "hisprof.h"
34
35#include "nobjmgr.h"
36#include "pistdimgapp.h"
37
38////////////////////////////////////////////////////////////////////
39////---------- Classe de fenetre de fit interactive ------------////
40////////////////////////////////////////////////////////////////////
41class PIAFitterWind : public PIWindow {
42public :
43 PIAFitterWind(PIStdImgApp* par, PIAFitter* fiter);
44 virtual ~PIAFitterWind();
45 virtual void Show();
46 virtual void Process(PIMessage msg, PIMsgHandler* sender, void* data=NULL);
47 void SetState(void);
48protected :
49 PIStdImgApp* mDap;
50 PIAFitter* mFitter;
51 PILabel** lab; int nlab;
52 PIText** txt; int ntxt;
53 PIButton** but; int nbut;
54 PICheckBox** ckb; int nckb;
55 PIOptMenu** pom; int npom;
56 bool ReFillGData;
57 int saveI1,saveI2, saveJ1, saveJ2;
58 double ErrValue;
59};
60
61////////////////////////////////////////////////////////////////////////
62//| --------------- Fit Lineaire a 1 et 2 dimensions ---------------
63//| Syntaxe:
64//| fitlin nom pnn [o:.aa,bb,cc, o;dd,ee,ff o:gg,hh,jj,kk etc...]
65//| avec:
66//| nom : cf commentaire ordre "fit"
67//| pnn : fit polynome degre nn avec classe Poly (lineaire) 1D ou 2D
68//| o:aa,...,bb : cf commentaires ordre "fit"
69//|
70////////////////////////////////////////////////////////////////////////
71//| --------------- Fit d'objets a 1 et 2 dimensions ---------------
72//| avec interface d'aide graphique
73//| Syntaxe:
74//| fitw nom func
75//| [p:p1,p2,...,pn s:s1,s2,...,sn m:m1,m2,...,mn M:M1,mM2,...,Mn f:f1,...,fn]
76//| [o:.aa,bb,cc, o;dd,ee,ff o:gg,hh,jj,kk etc...]
77//| cf commentaire ordre "fit"
78//|
79////////////////////////////////////////////////////////////////////////
80//| --------------- Fit d'objets a 1 et 2 dimensions ---------------
81//| Syntaxe:
82//| fit nom func
83//| [p:p1,p2,...,pn s:s1,s2,...,sn m:m1,m2,...,mn M:M1,mM2,...,Mn f:f1,...,fn]
84//| [o:.aa,bb,cc, o;dd,ee,ff o:gg,hh,jj,kk etc...]
85//|----- OBJET -----
86//| nom : nom de l'objet qui peut etre:
87//| fit-1D: Vector,Histo1D,HProf ou GeneraFitData(1D)
88//| fit-2D: Matrix,Histo2D,Image<T> ou GeneraFitData(2D)
89//|
90//|----- FUNCTION -----
91//| func : pnn : fit polynome degre nn avec GeneralFit (non-lineaire) 1D ou 2D
92//| : gnn : fit gaussienne (hauteur) + polynome de degre nn 1D
93//| : g : fit gaussienne (hauteur) 1D
94//| : enn : fit exponentielle + polynome de degre nn 1D
95//| : e : fit exponentielle 1D
96//| : Gnn : fit gaussienne (volume) + polynome de degre nn 1D
97//| : G : fit gaussienne (volume) 1D
98//| : : fit gaussienne+fond (volume) 2D
99//| : Gi : fit gaussienne+fond integree (volume) 2D
100//| : d : fit DL de gaussienne+fond (volume) 2D
101//| : di : fit DL de gaussienne+fond integree (volume) 2D
102//| : D : fit DL de gaussienne+fond avec coeff variable p6 (volume) 2D
103//| : Di : fit DL de gaussienne+fond integree avec coeff variable p6 (volume) 2D
104//| : M : fit Moffat+fond (expos=p6) (volume) 2D
105//| : Mi : fit Moffat+fond integree (expos=p6) (volume) 2D
106//| : Autre : fonction definie par l'utilisateur (cf commentaires ci apres)
107//|
108//|----- INIT PARAMETRES ET ETAT DU FIT -----
109//| p : p1,...,pn : valeur d'initialisation des parametres (def=0)
110//| s : s1,...,sn : valeur des steps de depart (def=1)
111//| m : m1,...,mn : valeur des minima (def=1)
112//| M : M1,...,Mn : valeur des maxima (def=-1) (max<=min : pas de limite)
113//| f : f1,...,fn : si >=1 parametre fixe sinon libre (def=0)
114//| - Remarque: si pi,si,mi ou Mi = '!' la valeur correspondante n'est pas changee
115//|
116//|----- OPTIONS -----
117//| o : options "o:Eaa.b,eaa.b,f,r,caa.b,Xaa.b"
118//| F : initialisation a partir des resultats et de l'etat du fit precedent
119//| (option non prioritaire sur les definitions p:,s:,m:,M:,f:o:)
120//| f : generation d'un Objet identique contenant la fonction fittee
121//| r : generation d'un Objet identique contenant les residus
122//| Xaa.b : aa.b valeur du DXi2 d'arret (def=1.e-3)
123//| Naa : aa nombre maximum d'iterations (def=100)
124//| la : niveau "a" de print: a=niveau de print Fit1/2D
125//| da : niveau "a" de debug: a=niveau de GeneralFit
126//| Ii1/i2 : numeros des bins X de l'histos utilises pour le fit [i1,i2]
127//|2D Jj1/j2 : numeros des bins Y de l'histos utilises pour le fit [j1,j2]
128//| - L'erreur est celle associee a l'objet (si elle existe),
129//| elle est mise a 1 sinon, sauf si E... ou e... est precise:
130//| Eaa.b : si |val|>=1 erreur = aa.b*sqrt(|val|)
131//| si |val|<1 erreur = aa.b
132//| si aa.b <=0 alors aa.b=1.0
133//| E seul est equivalent a E1.0
134//| eaa.b : erreur = aa.b
135//| si aa.b <=0 alors aa.b=1.0
136//| e seul est equivalent a e1.0
137//| xaa.b : demande de centrage: on fit x-aa.b au lieu de x)
138//| x : demande de centrage: on fit x-xc au lieu de x
139//| avec xc=abscisse du milieu de l'histogramme
140//| Actif pour exp+poly 1D, poly 1D
141//| Pour gauss+poly 1D, xc est le centre de la gaussienne.
142//|2D yaa.b et y : idem "xaa.b et x" mais pour y
143//|
144////////////////////////////////////////////////////////////////////////
145//| --------------- Def des fcts de fit par l'utilisateur ---------------
146//| Pour definir une fonction parametree utilisateur:
147//| Un exemple de fichier contenant des fonctions a fitter: PIext/userfitfunex.c
148//| > crfitfil nom_fichier nom_fonction nvar npar
149//| pour creer un squelette de fichier C contenant la fonction a fitter
150//| nom_fichier : nom du fichier ou est le code de la fonction
151//| c'est un fichier C de nom blabla.c
152//| nom_fonction : nom de la fonction ex: gauss2
153//| nvar : nombre de variables x[] (1D=1, 2D=2, etc...)
154//| npar : nombre de parametre p[]
155//| - La fonction est donc F(x[];p[]), soit par ex:
156//| double gauss2(double const* x,double const* p)
157//| - L'utilisateur a la possibilite de donner egalement la fonction
158//| retournant les derivees de la fonction par rapport aux parametres
159//| double gauss2_der(double const* x,double const* p,double* dp)
160//| ou dp[i] = dF(x[];p[])/dp[i]
161//| Par convention, cette fonction a le meme nom suivi de "_der"
162//| > crfitfun nom_fichier nom_fonction nvar npar
163//| Permet de compiler le fichier "nom_fichier" et de linker
164//| la fonction "nom_fonction"
165//| - Un meme fichier peut contenir plusieurs fonctions a fitter
166//|
167////////////////////////////////////////////////////////////////////////
168
169/* --Methode-- */
170PIAFitter::PIAFitter(PIACmd *piac, PIStdImgApp* app)
171: mApp(app), FWindFit(NULL)
172, mNObj(""), mObj(NULL), mGData(NULL)
173, mNPar(0), mNVar(0), mNBinX(0), mNBinY(0), mNData(0)
174, mPar(1), mStep(1), mMin(1), mMax(1), mFix(1)
175, mParSave(1), mStepSave(1), mMinSave(1), mMaxSave(1), mFixSave(1)
176, mFit(NULL)
177, mV(NULL), mH(NULL), mM(NULL), mH2(NULL), mIm(NULL), mG(NULL)
178, mFunction(NULL), mFName(""), mFunc(NULL)
179, mDlUFunc(NULL), mNameFitFunc(""), mFitFunc(NULL), mFitFuncDer(NULL)
180, mUFNVar(0), mUFNPar(0)
181{
182// Init et Set des options
183mPar=0.; mStep=1.; mMin=1.; mMax=-1.; mFix=0.;
184mParSave=0.; mStepSave=1.; mMinSave=1.; mMaxSave=-1.; mFixSave=0.;
185ResetOptions();
186mOptSave = mOpt;
187
188// enregistrement des ordres de fit
189string kw, usage;
190
191kw = "fit";
192usage = "Fitting function to DataObjects (Histo, Histo2D, Vector, ...)";
193usage += "\n Usage: fit nomobj func [Options]";
194usage += "\n [p:p1,...,pn s:s1,...,sn m:m1,...,mn M:M1,...,Mn o:... o:...]";
195piac->RegisterCommand(kw,usage,this,"Fitting");
196
197kw = "fitw";
198usage = "Fitting function to DataObjects with Interactive Window Help";
199usage += "\n Usage: fitw nomobj func [Options]";
200usage += "\n [p:p1,...,pn s:s1,...,sn m:m1,...,mn M:M1,...,Mn o:... o:...]";
201piac->RegisterCommand(kw,usage,this,"Fitting");
202
203kw = "fitlin";
204usage = "Linear Fitting of Polynoms to DataObjects";
205usage += "\n Usage: fitlin nomobj func [o:... o:...]";
206piac->RegisterCommand(kw,usage,this,"Fitting");
207
208kw = "crfitfun";
209usage = "Creation et link de function utilisateur pour le fit";
210usage += "\n Usage: crfitfun file func nvar npar";
211piac->RegisterCommand(kw,usage,this,"Fitting");
212
213kw = "crfitfil";
214usage = "Creation de fichier C pour function utilisateur du fit";
215usage += "\n Usage: crfitfil file func nvar npar";
216piac->RegisterCommand(kw,usage,this,"Fitting");
217}
218
219/* --Methode-- */
220PIAFitter::~PIAFitter()
221{
222if(mFit != NULL) {delete mFit; mFit=NULL;}
223if(mFunction != NULL) {delete mFunction; mFunction=NULL;}
224if(mFunc != NULL) {delete mFunc; mFunc=NULL;}
225if(mGData != NULL) {delete mGData; mGData=NULL;}
226if(FWindFit != NULL) {delete FWindFit; FWindFit=NULL;}
227if (mDlUFunc != NULL) {delete mDlUFunc; mDlUFunc=NULL;}
228}
229
230/* --Methode-- */
231void PIAFitter::ResetDPointer(void)
232// Reset des pointeurs sur les donnees
233{
234mNObj = ""; mObj=NULL;
235mV=NULL; mH=NULL; mM=NULL; mH2=NULL; mIm=NULL; mG=NULL;
236}
237
238/* --Methode-- */
239void PIAFitter::ResetOptions(void)
240// Reset des options
241{
242mOpt.okres = mOpt.okfun = false;
243mOpt.polcx = mOpt.polcy = 0;
244mOpt.xc = mOpt.yc = 0.;
245mOpt.stc2 = 1.e-3;
246mOpt.nstep = 100;
247mOpt.err_e = mOpt.err_E = -1.;
248mOpt.lp = 1; mOpt.lpg = 0;
249mOpt.i1 = mOpt.j1 = mOpt.i2 = mOpt.j2 = -1;
250mOpt.fromlastfit = false;
251}
252
253/* --Methode-- */
254void PIAFitter::DecodeOptions(string opt)
255// Decodage des options a partir des chaines de caracteres
256{
257ResetOptions();
258if(opt.length()<=0) return;
259
260size_t p,q; string dum;
261opt = "," + opt + ",";
262
263if(strstr(opt.c_str(),",F,")) // init options/param a partir du fit precedent
264 {mOpt = mOptSave; mOpt.fromlastfit = true;}
265if(strstr(opt.c_str(),",r,")) mOpt.okres = true; // residus
266if(strstr(opt.c_str(),",f,")) mOpt.okfun = true; // fonction fittee
267if(strstr(opt.c_str(),",x")) { // Demande de centrage (fit X=x-xc)
268 mOpt.polcx = 2; // Le centrage est calcule automatiquement
269 p = opt.find(",x"); q = opt.find_first_of(',',p+1);
270 dum = opt.substr(p,q-p);
271 if(dum.length()>2) {
272 sscanf(dum.c_str(),",x%lf",&mOpt.xc);
273 mOpt.polcx = 1; // Le centrage est fixe par la valeur lue
274 }
275}
276if(strstr(opt.c_str(),",y")) { // Demande de centrage (fit Y=y-yc)
277 mOpt.polcy = 2; // Le centrage est calcule automatiquement
278 p = opt.find(",y"); q = opt.find_first_of(',',p+1);
279 dum = opt.substr(p,q-p);
280 if(dum.length()>2) {
281 sscanf(dum.c_str(),",y%lf",&mOpt.yc);
282 mOpt.polcy = 1; // Le centrage est fixe par la valeur lue
283 }
284}
285if(strstr(opt.c_str(),",E")) { // Erreurs imposees a "sqrt(val)" ou "aa.b*sqrt(val)"
286 p = opt.find(",E"); q = opt.find_first_of(',',p+1);
287 dum = opt.substr(p,q-p);
288 if(dum.length()>2) sscanf(dum.c_str(),",E%lf",&mOpt.err_E);
289 if(mOpt.err_E<=0.) mOpt.err_E = 1.; //mOpt.err_e=-1.;
290}
291if(strstr(opt.c_str(),",e")) { // Erreurs imposees a "1" ou "aa.b"
292 p = opt.find(",e"); q = opt.find_first_of(',',p+1);
293 dum = opt.substr(p,q-p);
294 if(dum.length()>2) sscanf(dum.c_str(),",e%lf",&mOpt.err_e);
295 if(mOpt.err_e<=0.) mOpt.err_e = 1.; //mOpt.err_E=-1.;
296}
297if(strstr(opt.c_str(),",X")) { // Valeur du StopChi2
298 p = opt.find(",X"); q = opt.find_first_of(',',p+1);
299 dum = opt.substr(p,q-p);
300 if(dum.length()>2) sscanf(dum.c_str(),",X%lf",&mOpt.stc2);
301}
302if(strstr(opt.c_str(),",N")) { // Nombre maximum d'iterations
303 p = opt.find(",N"); q = opt.find_first_of(',',p+1);
304 dum = opt.substr(p,q-p);
305 if(dum.length()>2) sscanf(dum.c_str(),",N%d",&mOpt.nstep);
306}
307if(strstr(opt.c_str(),",l")) { // niveau de print
308 p = opt.find(",l"); q = opt.find_first_of(',',p+1);
309 dum = opt.substr(p,q-p);
310 if(dum.length()>2) sscanf(dum.c_str(),",l%d",&mOpt.lp);
311}
312if(strstr(opt.c_str(),",d")) { // niveau de debug du fit
313 p = opt.find(",d"); q = opt.find_first_of(',',p+1);
314 dum = opt.substr(p,q-p);
315 if(dum.length()>2) sscanf(dum.c_str(),",d%d",&mOpt.lpg);
316}
317if(strstr(opt.c_str(),",I")) { // intervalle de fit selon X
318 p = opt.find(",I"); q = opt.find_first_of(',',p+1);
319 dum = opt.substr(p,q-p);
320 if(dum.length()>2) sscanf(dum.c_str(),",I%d/%d",&mOpt.i1,&mOpt.i2);
321}
322if(strstr(opt.c_str(),",J")) { // intervalle de fit selon Y
323 p = opt.find(",J"); q = opt.find_first_of(',',p+1);
324 dum = opt.substr(p,q-p);
325 if(dum.length()>2) sscanf(dum.c_str(),",J%d/%d",&mOpt.j1,&mOpt.j2);
326}
327
328return;
329}
330
331/* --Methode-- */
332void PIAFitter::DecodeObject(string obj)
333// Decodage de l'objet a fitter
334{
335NamedObjMgr omg;
336AnyDataObj* mobj = omg.GetObj(obj);
337if(mobj == NULL)
338 {cout<<"PIAFitter::DecodeObject Error , Pas d'objet de nom "<<obj<<endl;
339 return;}
340
341// Reset
342ResetDPointer();
343mNVar = mNBinX = mNBinY = mNData = 0;
344
345mObj = mobj; mNObj = obj;
346string ctyp = typeid(*mObj).name();
347
348// 1D
349if (typeid(*mObj) == typeid(Vector)) {
350 mNVar = 1;
351 mV = (Vector*) mObj; mNBinX = mV->NElts(); mNBinY = 1;
352 }
353else if ( (typeid(*mObj) == typeid(HProf)) || (typeid(*mObj) == typeid(Histo)) ) {
354 mNVar = 1;
355 mH = (Histo*) mObj; mNBinX = mH->NBins(); mNBinY = 1;
356 }
357// 2D
358else if (typeid(*mObj) == typeid(Matrix)) {
359 mNVar = 2;
360 mM = (Matrix*) mObj; mNBinX = mM->NCol(); mNBinY =mM->NRows();
361 }
362else if (typeid(*mObj) == typeid(Histo2D)) {
363 mNVar = 2;
364 mH2 = (Histo2D*) mObj; mNBinX = mH2->NBinX(); mNBinY = mH2->NBinY();
365 }
366else if (typeid(*mObj) == typeid(GeneralFitData)) {
367 mG = (GeneralFitData*) mObj; mNBinX = mG->NData(); mNBinY = 1;
368 if( mG->NVar()==1) mNVar = 1;
369 else if(mG->NVar()==2) mNVar = 2;
370 else
371 {cout<<"PIAFitter::DecodeObject Error: \n"
372 <<" GeneralFitData ne peut avoir que 1 ou 2 variables d'abscisse: "
373 <<((GeneralFitData*) mObj)->NVar()<<endl;
374 return;}
375 }
376else if (dynamic_cast<RzImage*>(mObj)) {
377 mNVar = 2;
378 mIm = (RzImage*) mObj; mNBinX = mIm->XSize(); mNBinY = mIm->YSize();
379 }
380else {
381 cout<<"PIAFitter::DecodeObject Error , Objet n'est pas un "
382 <<"Histo1D/HProf/Vector/Histo2D/Image/Matrix/GeneralFitData "<<ctyp<<endl;
383 return;
384 }
385
386mNData = mNBinX*mNBinY;
387if(mNData<=0) {
388 cout<<"L'objet a "<<mNBinX<<","<<mNBinY<<" bins ("<<mNData<<")"<<endl;
389 return;
390}
391
392return;
393}
394
395/* --Methode-- */
396void PIAFitter::CheckOptions(void)
397// Check consistence des choix des options
398// ATTENTION: cette methode a besoin que l'objet a fitter soit deja lu
399{
400if(mOpt.lp<0) mOpt.lp = 0;
401if(mOpt.lpg<0) mOpt.lpg = 0;
402mOpt.i1 = (mOpt.i1<0||mOpt.i1>=mNBinX)? 0: mOpt.i1;
403mOpt.i2 = (mOpt.i2<0||mOpt.i2>=mNBinX||mOpt.i2<=mOpt.i1)? mNBinX-1: mOpt.i2;
404if(mNVar>=2) {
405 mOpt.j1 = (mOpt.j1<0||mOpt.j1>=mNBinY)? 0: mOpt.j1;
406 mOpt.j2 = (mOpt.j2<0||mOpt.j2>=mNBinY||mOpt.j2<=mOpt.j1)? mNBinY-1: mOpt.j2;
407} else mOpt.j2 = mOpt.j1 = 0;
408if(mOpt.polcx==2) {
409 if(mV||mM) mOpt.xc = (mOpt.i2-mOpt.i1+1)/2.;
410 else if(mH) mOpt.xc = (mH->XMin()+mH->XMax())/2.;
411 else if(mH2) mOpt.xc = (mH2->XMin()+mH2->XMax())/2.;
412 else if(mG) {double mini,maxi; mG->GetMinMax(2,mini,maxi); mOpt.xc=(mini+maxi)/2.;}
413 else if(mIm) {mOpt.xc = mIm->XOrg() * mIm->XPxSize()*(mOpt.i2-mOpt.i1+1)/2.;}
414}
415if(mOpt.polcy==2 && mNVar>=2) {
416 if(mM) mOpt.yc = (mOpt.j2-mOpt.j1+1)/2.;
417 if(mH2) mOpt.yc = (mH2->YMin()+mH2->YMax())/2.;
418 if(mG) {double mini,maxi; mG->GetMinMax(12,mini,maxi); mOpt.yc=(mini+maxi)/2.;}
419 if(mIm) {mOpt.yc = mIm->YOrg() * mIm->YPxSize()*(mOpt.j2-mOpt.j1+1)/2.;}
420}
421if(mOpt.err_e>0.) mOpt.err_E=-1.;
422if(mOpt.err_E>0.) mOpt.err_e =-1.;
423if(mOpt.stc2<=0.) mOpt.stc2 = 1.e-3;
424if(mOpt.nstep<2) mOpt.nstep = 100;
425
426return;
427}
428
429/* --Methode-- */
430void PIAFitter::PrintOptions(void)
431// Print des options
432{
433cout<<"Fit["<<mNBinX<<","<<mNBinY<<"] ("<<mNData<<") dim="<<mNVar<<":"
434 <<" Int=["<<mOpt.i1<<","<<mOpt.i2<<"],["<<mOpt.j1<<","<<mOpt.j2<<"]"<<endl
435 <<" Cent="<<mOpt.polcx<<","<<mOpt.polcy<<","<<mOpt.xc<<"+x"<<","<<mOpt.yc<<"+y"
436 <<" TypE="<<mOpt.err_e<<","<<mOpt.err_E
437 <<" StpX2="<<mOpt.stc2<<" Nstep="<<mOpt.nstep<<endl
438 <<" Init.LFit="<<mOpt.fromlastfit
439 <<" lp,lpg="<<mOpt.lp<<","<<mOpt.lpg<<endl;
440return;
441}
442
443/* --Methode-- */
444void PIAFitter::FillGData(void)
445// Fill generalfitdata pour le fit
446// ATTENTION: cette methode a besoin que les options soient decodees et checkees
447{
448if(mNData<=0) return;
449if(mGData!=NULL) {delete mGData; mGData=NULL;}
450mGData = new GeneralFitData(mNVar,mNData,0);
451
452for(int i=mOpt.i1;i<=mOpt.i2;i++) for(int j=mOpt.j1;j<=mOpt.j2;j++) {
453 double x,y=0.,f,e;
454
455 if(mV)
456 {x= (double) i; f=(*mV)(i); e=1.;}
457 else if(mH)
458 {x=mH->BinCenter(i); f=(*mH)(i); e=(mH->HasErrors())?mH->Error(i):1.;}
459 else if(mM)
460 {x=(double) i; y=(double) j; f=(*mM)(j,i); e=1.;}
461 else if(mH2)
462 {float xf,yf; mH2->BinCenter(i,j,xf,yf); x=(double)xf; y=(double)yf;
463 f=(*mH2)(i,j); e=(mH2->HasErrors())?mH2->Error(i,j):1.;}
464 else if(mIm)
465 {x=mIm->XOrg()+(i+0.5)*mIm->XPxSize(); y=mIm->YOrg()+(j+0.5)*mIm->YPxSize();
466 f=mIm->DValue(i,j); e=1.;}
467 else if(mG&&mNVar==1) {x= mG->X(i); f=mG->Val(i); e=mG->EVal(i);}
468 else if(mG&&mNVar==2) {x= mG->X(i); y=mG->Y(i); f=mG->Val(i); e=mG->EVal(i);}
469 else x=y=f=e=0.;
470
471 // Gestion des erreurs a utiliser
472 if(mOpt.err_e>0.) e=mOpt.err_e;
473 else if(mOpt.err_E>0.) {e=(f<-1.||f>1.)?mOpt.err_E*sqrt(fabs(f)):mOpt.err_E;}
474
475 // Remplissage de generalfitdata
476 if(mFName[0]=='p') {x -= mOpt.xc; if(mNVar>=2) y -= mOpt.yc;}
477 if(mNVar==1) mGData->AddData1(x,f,e);
478 else if(mNVar==2) mGData->AddData2(x,y,f,e);
479}
480
481if(mGData->NData()<=0) {
482 cout<<"Pas de donnees dans GeneralFitData: "<<mGData->NData()<<endl;
483 return;
484}
485
486if(mOpt.lpg>1) {
487 mGData->PrintStatus();
488 mGData->PrintData(0);
489 mGData->PrintData(mGData->NData()-1);
490}
491
492return;
493}
494
495/* --Methode-- */
496void PIAFitter::DecodeFunction(string func)
497// Fonction a fitter
498// ATTENTION: cette methode a besoin que les donnees soient lues
499{
500if(func.length()<=0)
501 {cout<<"PIAFitter::DecodeFunction Donnez un nom de fonction a fitter."<<endl;
502 return;}
503if(mFunction!=NULL) {delete mFunction; mFunction=NULL;}
504if(mFunc!=NULL) {delete mFunc; mFunc=NULL;}
505mNPar=0;
506
507if(func == mNameFitFunc) { // Fonction utilisateur
508 mFunc = new GeneralFunc(mUFNVar,mUFNPar,mFitFunc,mFitFuncDer);
509 mNPar = mUFNPar;
510 cout<<"Fonction utilisateur nvar="<<mUFNVar<<", npar="<<mNPar<<endl;
511} else if(func[0]=='p' && mNVar==1) { //polynome
512 int degre=0; if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
513 cout<<"Fit polynome 1D de degre "<<degre<<endl;
514 Polyn1D* myf = new Polyn1D(degre,mOpt.xc);
515 mFunction = myf; mNPar = mFunction->NPar(); mFName = func;
516
517} else if(func[0]=='e' && mNVar==1) { //exponentielle
518 int degre=-1; if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
519 cout<<"Fit d'exponentielle+polynome 1D de degre "<<degre<<endl;
520 Exp1DPol* myf;
521 if(degre>=0) myf = new Exp1DPol((unsigned int)degre,mOpt.xc);
522 else myf = new Exp1DPol(mOpt.xc);
523 mFunction = myf; mNPar = mFunction->NPar(); mFName = func;
524
525} else if(func[0]=='g' && mNVar==1) { //gaussienne en hauteur
526 int degre=-1; if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
527 cout<<"Fit de Gaussienne_en_hauteur+polynome 1D de degre "<<degre<<endl;
528 Gauss1DPol* myf;
529 if(degre>=0) myf = new Gauss1DPol((unsigned int)degre,((mOpt.polcx)?true:false));
530 else { bool bfg = (mOpt.polcx)?true:false; myf = new Gauss1DPol(bfg); }
531 mFunction = myf; mNPar = mFunction->NPar(); mFName = func;
532
533} else if(func[0]=='G' && mNVar==1) { //gaussienne en volume
534 int degre=-1; if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
535 cout<<"Fit de Gaussienne_en_volume+polynome 1D de degre "<<degre<<endl;
536 GaussN1DPol* myf;
537 if(degre>=0) myf = new GaussN1DPol((unsigned int)degre,((mOpt.polcx)?true:false));
538 else { bool bfg = (mOpt.polcx)?true:false; myf = new GaussN1DPol(bfg); }
539 mFunction = myf; mNPar = mFunction->NPar(); mFName = func;
540
541} else if(func[0]=='p' && mNVar==2) { //polynome 2D
542 int degre=0; if(func.length()>1) sscanf(func.c_str()+1,"%d",&degre);
543 cout<<"Fit polynome 2D de degre "<<degre<<endl;
544 Polyn2D* myf = new Polyn2D(degre,mOpt.xc,mOpt.yc);
545 mFunction = myf; mNPar = mFunction->NPar(); mFName = func;
546
547} else if(func[0]=='G' && mNVar==2) { //gaussienne+fond en volume
548 int integ=0; if(func.length()>1) if(func[1]=='i') integ=1;
549 cout<<"Fit de Gaussienne+Fond 2D integ="<<integ<<endl;
550 if(integ) {GauRhInt2D* myf = new GauRhInt2D; mFunction = myf;}
551 else {GauRho2D* myf = new GauRho2D; mFunction = myf;}
552 mNPar = mFunction->NPar(); mFName = func;
553
554} else if(func[0]=='d' && mNVar==2) { //DL gaussienne+fond en volume
555 int integ=0; if(func.length()>1) if(func[1]=='i') integ=1;
556 cout<<"Fit de DL de Gaussienne+Fond 2D integ="<<integ<<endl;
557 if(integ) {GdlRhInt2D* myf = new GdlRhInt2D; mFunction = myf;}
558 else {GdlRho2D* myf = new GdlRho2D; mFunction = myf;}
559 mNPar = mFunction->NPar(); mFName = func;
560
561} else if(func[0]=='D' && mNVar==2) { //DL gaussienne+fond avec coeff variable p6 en volume
562 int integ=0; if(func.length()>1) if(func[1]=='i') integ=1;
563 cout<<"Fit de DL de Gaussienne+Fond avec coeff variable (p6) 2D integ="<<integ<<endl;
564 if(integ) {Gdl1RhInt2D* myf = new Gdl1RhInt2D; mFunction = myf;}
565 else {Gdl1Rho2D* myf = new Gdl1Rho2D; mFunction = myf;}
566 mNPar = mFunction->NPar(); mFName = func;
567
568} else if(func[0]=='M' && mNVar==2) { //Moffat+fond (volume)
569 int integ=0; if(func.length()>1) if(func[1]=='i') integ=1;
570 cout<<"Fit de Moffat+Fond (expos=p6) 2D integ="<<integ<<endl;
571 if(integ) {MofRhInt2D* myf = new MofRhInt2D; mFunction = myf;}
572 else {MofRho2D* myf = new MofRho2D; mFunction = myf;}
573 mNPar = mFunction->NPar(); mFName = func;
574
575} else {
576 cout<<"Fonction "<<func<<" inconnue pour la dim "<<mNVar<<endl;
577 return;
578}
579
580return;
581}
582
583/* --Methode-- */
584void PIAFitter::ReSetParam(void)
585// Pour (re)dimensionner les tableaux des parametres
586// ATTENTION: cette methode a besoin que la fonction soit connue
587{
588if(mNPar==0) return;
589mPar.Realloc(mNPar); mPar=0.;
590mStep.Realloc(mNPar); mStep=1.;
591mMin.Realloc(mNPar); mMin=1.;
592mMax.Realloc(mNPar); mMax=-1.;
593mFix.Realloc(mNPar); mFix=0.;
594return;
595}
596
597/* --Methode-- */
598void PIAFitter::InitParFromLastFit(void)
599// Initialisation du fit a partir des resultats du fit precedent
600// (Les arguments de la ligne de commande sont prioritaires)
601// ATTENTION: cette methode a besoin que la fonction soit connue
602{
603if(mNPar==0) return;
604int n = mParSave.NElts();
605if(n==0) return;
606
607for(int i=0;i<mNPar;i++) {
608 if(i>=n) break;
609 mPar(i) = mParSave(i);
610 mStep(i) = mStepSave(i);
611 mMin(i) = mMinSave(i);
612 mMax(i) = mMaxSave(i);
613 mFix(i) = mFixSave(i);
614}
615
616return;
617}
618
619/* --Methode-- */
620void PIAFitter::DecodeParam(string par,string step,string min,string max,string fix)
621// Decodage des parametres
622// ATTENTION: cette methode a besoin que la fonction soit connue
623// et que les options soient decodees
624{
625if(mNPar==0) return;
626ReSetParam();
627if(mOpt.fromlastfit) InitParFromLastFit();
628Vector* v=NULL; string* s=NULL;
629for(int j=0;j<5;j++) {
630 if(j==0) {v=&mPar; s=&par;} else if(j==1) {v=&mStep; s=&step;}
631 else if(j==2) {v=&mMin; s=&min;} else if(j==3) {v=&mMax; s=&max;}
632 else if(j==4) {v=&mFix; s=&fix;}
633 if(s->length()>0) *s += ",";
634 for(int i=0;i<mNPar;i++) {
635 if(s->length()<=0) break;
636 if((*s)[0]!='!') sscanf(s->c_str(),"%lf",&(*v)(i));
637 size_t p = s->find_first_of(',') + 1;
638 if(p>=s->length()) *s = ""; else *s = s->substr(p);
639 }
640}
641return;
642}
643
644/* --Methode-- */
645int PIAFitter::DoFit(void)
646// Fit avec generalfit
647{
648if(mNPar==0) return -1000;
649if((mFunction==NULL && mFunc==NULL) || mGData==NULL) {
650 cout<<"PIAFitter::DoFit error: mFunction,mFunc="<<mFunction<<","<<mFunc
651 <<" mGData="<<mGData<<endl;
652 return -1001;
653}
654// classe GeneraFit
655if(mFit != NULL) {delete mFit; mFit= NULL;}
656if(mFunction) mFit= new GeneralFit(mFunction);
657 else mFit= new GeneralFit(mFunc);
658// Options et Set de GeneralFit
659mFit->SetDebug(mOpt.lpg);
660mFit->SetData(mGData);
661mFit->SetStopChi2(mOpt.stc2);
662mFit->SetMaxStep(mOpt.nstep);
663{for(int i=0;i<mNPar;i++) {
664 char str[10];
665 sprintf(str,"P%d",i);
666 mFit->SetParam(i,str,mPar(i),mStep(i),mMin(i),mMax(i));
667 if((int)(mFix(i)+0.001)>=1.) mFit->SetFix(i);
668}}
669if(mOpt.lp>1) mFit->PrintFit();
670
671// Fit
672mParSave = mPar;
673mStepSave = mStep; mMinSave = mMin; mMaxSave = mMax; mFixSave = mFix;
674mOptSave = mOpt;
675double c2r = -1.;
676int rcfit = mFit->Fit();
677if(mOpt.lp>0) mFit->PrintFit();
678if(rcfit>0) {
679 c2r = mFit->GetChi2Red();
680 cout<<"C2r_Reduit = "<<c2r<<" nstep="<<mFit->GetNStep()<<" rc="<<rcfit<<endl;
681 for(int i=0;i<mNPar;i++) mParSave(i)=mFit->GetParm(i);
682} else {
683 cout<<"echec Fit, rc = "<<rcfit<<" nstep="<<mFit->GetNStep()<<endl;
684 mFit->PrintFitErr(rcfit);
685}
686
687return rcfit;
688}
689
690/* --Methode-- */
691void PIAFitter::FitFunRes(void)
692// Mise a disposition des resultats (fonction fitee et residus)
693{
694if(mFit && (mOpt.okres||mOpt.okfun)) {
695 NamedObjMgr omg;
696 string nomres = mNObj + "_res";
697 string nomfun = mNObj + "_fun";
698 if(mV) {
699 if(mOpt.okres) {Vector* ob = mV->FitResidus(*mFit); if(ob) omg.AddObj(ob,nomres);}
700 if(mOpt.okfun) {Vector* ob = mV->FitFunction(*mFit); if(ob) omg.AddObj(ob,nomfun);}
701 } else if(mH) {
702 if(mOpt.okres) {Histo* ob = mH->FitResidus(*mFit); if(ob) omg.AddObj(ob,nomres);}
703 if(mOpt.okfun) {Histo* ob = mH->FitFunction(*mFit); if(ob) omg.AddObj(ob,nomfun);}
704 } else if(mM) {
705 if(mOpt.okres) {Matrix* ob = mM->FitResidus(*mFit); if(ob) omg.AddObj(ob,nomres);}
706 if(mOpt.okfun) {Matrix* ob = mM->FitFunction(*mFit); if(ob) omg.AddObj(ob,nomfun);}
707 } else if(mH2) {
708 if(mOpt.okres) {Histo2D* ob = mH2->FitResidus(*mFit); if(ob) omg.AddObj(ob,nomres);}
709 if(mOpt.okfun) {Histo2D* ob = mH2->FitFunction(*mFit); if(ob) omg.AddObj(ob,nomfun);}
710 } else if(mIm) {
711 if(mOpt.okres) {
712 RzImage* ob = mIm->FitResidus(*mFit);
713 if(ob) {ImageR4* rob = new ImageR4(*ob); omg.AddObj(rob,nomres); delete ob;}
714 }
715 if(mOpt.okfun) {
716 RzImage* ob = mIm->FitFunction(*mFit);
717 if(ob) {ImageR4* rob = new ImageR4(*ob); omg.AddObj(rob,nomfun); delete ob;}
718 }
719 } else if(mG) {
720 if(mOpt.okres) {GeneralFitData* ob = mG->FitResidus(*mFit); if(ob) omg.AddObj(ob,nomres);}
721 if(mOpt.okfun) {GeneralFitData* ob = mG->FitFunction(*mFit); if(ob) omg.AddObj(ob,nomfun);}
722 }
723}
724return;
725}
726
727/* --Methode-- */
728void PIAFitter::LinFit(void)
729// Fit lineaire de polynomes
730{
731if(mFName.length()<=0)
732 {cout<<"PIAFitter::LinFit Donnez un nom de fonction a fitter."<<endl;
733 return;}
734if(!mGData)
735 {cout<<"PIAFitter::LinFit pas de Data."<<endl; return;}
736if(mFName[0]!='p')
737 {cout<<"PIAFitter::LinFit fonction non prevue "<<mFName<<endl;
738 return;}
739int degre = 0;
740if(mFName.length()>1) sscanf(mFName.c_str()+1,"%d",&degre);
741
742if(mNVar==1) { // Fit lineaire de polynome 1D
743 cout<<"Fit (lineaire) 1D polynome de degre "<<degre<<endl;
744 Poly p1(0);
745 double c2rl = mGData->PolFit(0,p1,degre);
746 cout<<"C2r_lineaire = "<<c2rl<<endl;
747 if(mOpt.lp>0) cout<<p1<<endl;
748} else if(mNVar==2) { // Fit lineaire de polynome 2D
749 cout<<"Fit (lineaire) polynome 2D de degre "<<degre<<endl;
750 Poly2 p2(0);
751 double c2rl = mGData->PolFit(0,1,p2,degre);
752 cout<<"C2r_lineaire = "<<c2rl<<endl;
753 if(mOpt.lp>0) cout<<p2<<endl;
754}
755
756return;
757}
758
759/* --Methode-- */
760void PIAFitter::LinkFitFunc(string const& fname, string const& funcname
761 ,int nvar,int npar)
762// Compilation et link dynamique de la fonction a fitter
763{
764if(nvar<1 || npar<1) return;
765if(mDlUFunc != NULL) {delete mDlUFunc; mDlUFunc=NULL;}
766mFitFunc = NULL; mFitFuncDer = NULL;
767mUFNVar= mUFNPar = 0; mNameFitFunc = "";
768
769// On compile le fichier C, on fabrique un .so et on l'ouvre
770mDlUFunc = PDynLinkMgr::BuildFromCFile(fname);
771if(mDlUFunc == NULL)
772 {cerr<<"PIAFitter::LinkFitFunc() Erreur building .so "<<endl; return;}
773
774// On recupere la fonction de fit
775mFitFunc = (DlUserFitFunc) mDlUFunc->GetFunction(funcname);
776if(mFitFunc == NULL) // Probleme, on n'a pas trouve etc ...
777 {cerr<<"PIAFitter::LinkFitFunc() Erreur linking function"<<funcname<<endl;
778 delete mDlUFunc; mDlUFunc = NULL; return;}
779
780// On recupere la fonction des derivees
781string fnameder = funcname+"_der";
782mFitFuncDer = (DlUserFitFuncDer) mDlUFunc->GetFunction(fnameder);
783if(mFitFuncDer == NULL)
784 {cout<<"PIAFitter::LinkFitFunc() No derive for"<<fnameder<<endl;}
785
786// Variables et pointeurs definissant le link a la fonction utilisateur
787mNameFitFunc = funcname;
788mUFNVar = nvar;
789mUFNPar = npar;
790cout<<"PIAFitter::LinkFitFunc() UserFitFunc="<<mNameFitFunc
791 <<" NVar="<<mUFNVar<<" NPar="<<mUFNPar<<endl;
792
793return;
794}
795
796/* --Methode-- */
797void PIAFitter::PrepareFitFuncCFile(string const& fname, string const& funcname
798 , int nvar, int npar)
799// Pour ecrire un squelette des fichier C avec pour une fonction a fitter
800{
801if(fname.length()<=0 || funcname.length()<=0 || nvar<=0 || npar <=0) {
802 cerr<<"PIAFitter::PrepareFitFuncCFile: error "
803 <<" nvar="<<nvar<<" npar="<<npar<<endl
804 <<" file_name="<<fname<<" func_name="<<funcname<<endl
805 <<" ---- No File Created ! ----"<<endl;
806 return;
807}
808
809FILE* fip = fopen(fname.c_str(),"w");
810if(fip==NULL)
811 {cerr<<"PIAFitter::PrepareFitFuncCFile: error opening file for writing\n"
812 <<fname<<endl; return;}
813
814fprintf(fip,"#include <stdio.h>\n#include <stdlib.h>\n#include <math.h>\n\n");
815fprintf(fip,"/* - many functions can be written in this file: */\n");
816fprintf(fip,"/* here is an exemple for creation of function %s */\n\n"
817 ,funcname.c_str());
818fprintf(fip,"double %s(double const* x,double const* p)\n",funcname.c_str());
819fprintf(fip,"/* Function : F(x[0-%d];p[0-%d]) */\n",nvar-1,npar-1);
820fprintf(fip,"/* variables x1=x[0],x2=x[1],...,x%d=x[%d] */\n",nvar,nvar-1);
821fprintf(fip,"/* parameters p1=p[0],p2=p[1],...,p%d=p[%d] */\n",nvar,nvar-1);
822fprintf(fip,"{\n val = F(x[];p[]);\n return(val);\n}\n\n");
823string fder = funcname + "_der";
824fprintf(fip,"/* That function is optional */\n");
825fprintf(fip,"/*\n");
826fprintf(fip,"double %s(double const* x,double const* p,double* dp)\n",fder.c_str());
827fprintf(fip,"{\n");
828for(int i=0;i<npar;i++) fprintf(fip," dp[%d] = dF(x[];p[])/dp[%d]\n",i,i);
829fprintf(fip," val = F(x[];p[]);\n return(val);\n}\n\n");
830fprintf(fip,"*/\n");
831
832fclose(fip);
833}
834
835/* --Methode-- */
836int PIAFitter::Execute(string& kw, vector<string>& tokens)
837// Execution des commandes de fit : "fit" et "fitw"
838{
839////// Les commandes pour creer les fonctions de fit utilisateur
840if(kw == "crfitfil") {
841 if (tokens.size() < 4)
842 {cout<<"Usage:crfitfil filename funcname nvar npar"<<endl;return(-1);}
843 int nvar = atoi(tokens[2].c_str()), npar = atoi(tokens[3].c_str());
844 PrepareFitFuncCFile(tokens[0],tokens[1],nvar,npar);
845 return(0);
846} else if (kw == "crfitfun") {
847 if (tokens.size() < 4)
848 {cout<<"Usage:crfitfun filename funcname nvar npar"<<endl;return(-1);}
849 int nvar = atoi(tokens[2].c_str()), npar = atoi(tokens[3].c_str());
850 LinkFitFunc(tokens[0],tokens[1],nvar,npar);
851 return(0);
852}
853
854////// Les commandes pour fitter
855if(tokens.size() < 2) {
856 cout
857 <<"Usage:fit nomobj func [Options]\n"
858 <<" ou fitw nomobj func [Options]\n"
859 <<" avec Options: [p:p1,...,pn s:s1,...,sn m:m1,...,mn \n"
860 <<" M:M1,...,Mn f:f1,...,fn o:.,.,. o:.,.,.]\n"
861 <<" ou fitlin nomobj func [o:.,.,. o:.,.,.]\n";
862 return(-1);
863}
864// decodage des arguments de la commande
865string mSp,mSs,mSm,mSM,mSf,mSo;
866mSp=""; mSs=""; mSm=""; mSM=""; mSf=""; mSo="";
867if(tokens.size()>2)
868 for(int ip=2;ip<(int)tokens.size();ip++) {
869 if(tokens[ip].length()<=2) continue;
870 const char *c = tokens[ip].c_str();
871 if(c[1]!=':') continue;
872 if(c[0]=='p') mSp=c+2;
873 else if(c[0]=='s') mSs=c+2;
874 else if(c[0]=='m') mSm=c+2;
875 else if(c[0]=='M') mSM=c+2;
876 else if(c[0]=='f') mSf=c+2;
877 else if(c[0]=='o') {mSo += ","; mSo += c+2;}
878 }
879// Execution des commandes
880DecodeOptions(mSo);
881DecodeObject(tokens[0]);
882CheckOptions();
883if(mOpt.lp>0) PrintOptions();
884FillGData();
885DecodeFunction(tokens[1]);
886DecodeParam(mSp,mSs,mSm,mSM,mSf);
887if (kw == "fit") {
888 int rc = DoFit();
889 if(rc>=0) FitFunRes();
890} else if (kw == "fitlin") {
891 LinFit();
892} else if (kw == "fitw") {
893 if(FWindFit!=NULL) {delete FWindFit; FWindFit=NULL;}
894 FWindFit = new PIAFitterWind(mApp,this);
895 FWindFit->Show();
896}
897
898return(0);
899}
900
901////////////////////////////////////////////////////////////////////
902////---------- Classe de fenetre de fit interactive ------------////
903////////////////////////////////////////////////////////////////////
904
905/* --Methode-- */
906PIAFitterWind::PIAFitterWind(PIStdImgApp* par, PIAFitter* fiter)
907: PIWindow((PIMsgHandler*)par, "PIAFitter", PIWK_normal, 240, 240, 150, 150)
908, mDap(par), mFitter(fiter), ReFillGData(false)
909, saveI1(0), saveI2(0), saveJ1(0), saveJ2(0), ErrValue(1.)
910{
911COMPATIBLE_ASSERT(mFitter);
912int npar = mFitter->mNPar;
913
914// Alloc de la taille
915nlab = 14+npar; lab = new PILabel*[nlab]; {for(int i=0;i<nlab;i++) lab[i]=NULL;}
916ntxt = 10+4*npar; txt = new PIText*[ntxt]; {for(int i=0;i<ntxt;i++) txt[i]=NULL;}
917nbut = 5; but = new PIButton*[nbut]; {for(int i=0;i<nbut;i++) but[i]=NULL;}
918nckb = 2+npar; ckb = new PICheckBox*[nckb]; {for(int i=0;i<nckb;i++) ckb[i]=NULL;}
919npom = 3; pom = new PIOptMenu*[npom]; {for(int i=0;i<npom;i++) pom[i]=NULL;}
920
921// On definit la taille a partir de la taille par defaut des composantes
922// bsx,bsy = environ 6 lettres
923int bsx,bsy;
924par->PrefCompSz(bsx, bsy);
925int spx = (bsx>=10) ? bsx/10 : 1; // intervalle entre lettres X
926int spy = (bsy>=6) ? bsy/6 : 1; // intervalle entre lettres Y
927
928int wszx = (5*bsx+5*spx)+spx;
929int wszy = 8*(bsy+spy)+npar*(bsy+spy)+2*spy;
930SetSize(wszx, wszy);
931int cpx,cpy;
932
933// new ligne
934cpx = spx; cpy = spy;
935string dum = "Object: "+mFitter->mNObj+" , Fun: "+mFitter->mFName;
936lab[0] = new PILabel(this,dum.c_str(),4*bsx,bsy,cpx,cpy);
937lab[0]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
938
939// new ligne
940cpx = spx; cpy += bsy+spy;
941but[0] = new PIButton(this,"Fr.Last",111,bsx,bsy,cpx,cpy);
942but[0]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
943cpx += bsx+spx;
944but[1] = new PIButton(this,"Default",222,bsx,bsy,cpx,cpy);
945but[1]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
946cpx += bsx+spx;
947ckb[0] = new PICheckBox(this,"Gen Fun",1001,bsx,bsy,cpx,cpy);
948ckb[0]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
949cpx += bsx+2*spx;
950ckb[1] = new PICheckBox(this,"Res",1002,bsx,bsy,cpx,cpy);
951ckb[1]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
952
953// new ligne
954cpx = spx; cpy += bsy+spy;
955pom[0] = new PIOptMenu(this,"Err def",1.5*bsx,bsy,cpx,cpy);
956 pom[0]->AppendItem("Err def",2101);
957 pom[0]->AppendItem("Err cste",2102);
958 pom[0]->AppendItem("Err sqrt",2103);
959 pom[0]->SetValue(2101);
960pom[0]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
961cpx += int(1.5*bsx+spx);
962txt[0] = new PIText(this,"",bsx,bsy,cpx,cpy);
963txt[0]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
964cpx += bsx+spx;
965lab[1] = new PILabel(this,"Prt,dbg",bsx,bsy,cpx,cpy);
966lab[1]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
967cpx += bsx+spx;
968txt[1] = new PIText(this,"",bsx,bsy,cpx,cpy);
969txt[1]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
970
971// new ligne
972cpx = spx; cpy += bsy+spy;
973lab[2] = new PILabel(this,"Stop Xi2",1.25*bsx,bsy,cpx,cpy);
974lab[2]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
975cpx += int(1.25*bsx+spx);
976txt[2] = new PIText(this,"",bsx,bsy,cpx,cpy);
977txt[2]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
978cpx += bsx+spx;
979lab[3] = new PILabel(this,"Iter",bsx, bsy, cpx, cpy);
980lab[3]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
981cpx += bsx+spx;
982txt[3] = new PIText(this,"",bsx,bsy,cpx,cpy);
983txt[3]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
984
985// new ligne
986cpx = spx; cpy += bsy+spy;
987lab[4] = new PILabel(this,"Range X",bsx,bsy,cpx,cpy);
988lab[4]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
989cpx += bsx+spx;
990txt[4] = new PIText(this,"",0.75*bsx,bsy,cpx,cpy);
991txt[4]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
992cpx += int(0.75*bsx+spx);
993txt[5] = new PIText(this,"",0.75*bsx,bsy,cpx,cpy);
994txt[5]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
995cpx += int(0.75*bsx+spx);
996lab[5] = new PILabel(this,"Y",bsx/2,bsy,cpx,cpy);
997lab[5]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
998cpx += int(bsx/2.+spx);
999txt[6] = new PIText(this,"",0.75*bsx,bsy,cpx,cpy);
1000txt[6]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1001cpx += int(0.75*bsx+spx);
1002txt[7] = new PIText(this,"",0.75*bsx,bsy,cpx,cpy);
1003txt[7]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1004
1005// new ligne
1006cpx = spx; cpy += bsy+spy;
1007lab[6] = new PILabel(this,"Cent",bsx/2,bsy,cpx,cpy);
1008lab[6]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1009cpx += int(bsx/2.+spx);
1010pom[1] = new PIOptMenu(this,"X.No",1.25*bsx,bsy,cpx,cpy);
1011 pom[1]->AppendItem("X.No",2001);
1012 pom[1]->AppendItem("X.Value",2002);
1013 pom[1]->AppendItem("X.Auto",2003);
1014 pom[1]->SetValue(2001);
1015pom[1]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1016cpx += int(1.25*bsx+spx);
1017txt[8] = new PIText(this,"",bsx,bsy,cpx,cpy);
1018txt[8]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1019cpx += bsx+spx;
1020pom[2] = new PIOptMenu(this,"Y.No",1.25*bsx,bsy,cpx,cpy);
1021 pom[2]->AppendItem("Y.No",2011);
1022 pom[2]->AppendItem("Y.Value",2012);
1023 pom[2]->AppendItem("Y.Auto",2013);
1024 pom[2]->SetValue(2011);
1025pom[2]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1026cpx += int(1.25*bsx+spx);
1027txt[9] = new PIText(this,"",bsx,bsy,cpx,cpy);
1028txt[9]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1029
1030// new ligne
1031cpx = spx; cpy += bsy+spy;
1032but[2] = new PIButton(this, "Dismiss",777,bsx,bsy,cpx,cpy);
1033but[2]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1034cpx += bsx+spx;
1035but[3] = new PIButton(this, "Check",666,bsx,bsy,cpx,cpy);
1036but[3]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1037cpx += bsx+spx;
1038but[4] = new PIButton(this, "Fit",333,bsx,bsy,cpx,cpy);
1039but[4]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1040cpx += bsx+spx;
1041lab[7] = new PILabel(this,mFitter->mFName.c_str(),2*bsx,bsy,cpx,cpy);
1042lab[7]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1043
1044// new ligne
1045cpx = spx; cpy += bsy+spy;
1046lab[8] = new PILabel(this,"Par",bsx/2,bsy,cpx,cpy);
1047lab[8]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1048cpx += int(bsx/2.+spx);
1049lab[9] = new PILabel(this,"Fx",bsx/3,bsy,cpx,cpy);
1050lab[9]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1051cpx += int(bsx/3.+spx);
1052lab[10] = new PILabel(this,"Init",bsx,bsy,cpx,cpy);
1053lab[10]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1054cpx += bsx+spx;
1055lab[11] = new PILabel(this,"Step",bsx,bsy,cpx,cpy);
1056lab[11]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1057cpx += bsx+spx;
1058lab[12] = new PILabel(this,"Min",bsx,bsy,cpx,cpy);
1059lab[12]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1060cpx += bsx+spx;
1061lab[13] = new PILabel(this,"Max",bsx,bsy,cpx,cpy);
1062lab[13]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1063
1064// new lines: Parametres (npar lignes)
1065if(npar>0) {
1066 for(int i=0;i<npar;i++) {
1067 cpx = spx; cpy += bsy+spy;
1068 char str[8]; sprintf(str,"P%d",i);
1069 lab[14+i] = new PILabel(this,str,bsx/2,bsy,cpx,cpy);
1070 lab[14+i]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1071 cpx += int(bsx/2.+spx);
1072 // Attention pas de message entre 1500 et 2000 (permet 500 parametres!)
1073 ckb[2+i] = new PICheckBox(this,"",1500+i,bsx/3,bsy,cpx,cpy+bsy/2);
1074 ckb[2+i]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1075 cpx+=int(bsx/3.+spx);
1076 for(int j=0;j<4;j++) {
1077 txt[10+4*i+j] = new PIText(this,"",bsx,bsy,cpx,cpy);
1078 txt[10+4*i+j]->SetBinding(PIBK_elastic,PIBK_elastic, PIBK_elastic,PIBK_elastic);
1079 cpx += bsx+spx;
1080 }
1081 }
1082}
1083
1084return;
1085}
1086
1087/* --Methode-- */
1088PIAFitterWind::~PIAFitterWind()
1089{
1090int i;
1091for(i=0;i<nlab;i++) if(lab[i]) {delete lab[i]; lab[i]=NULL;}
1092for(i=0;i<ntxt;i++) if(txt[i]) {delete txt[i]; txt[i]=NULL;}
1093for(i=0;i<nbut;i++) if(but[i]) {delete but[i]; but[i]=NULL;}
1094for(i=0;i<nckb;i++) if(ckb[i]) {delete ckb[i]; ckb[i]=NULL;}
1095for(i=0;i<npom;i++) if(pom[i]) {delete pom[i]; pom[i]=NULL;}
1096if(lab) {delete lab; lab=NULL; nlab=0;}
1097if(txt) {delete txt; txt=NULL; ntxt=0;}
1098if(but) {delete but; but=NULL; nbut=0;}
1099if(ckb) {delete ckb; ckb=NULL; nckb=0;}
1100if(pom) {delete pom; pom=NULL; npom=0;}
1101}
1102
1103/* --Methode-- */
1104void PIAFitterWind::Show()
1105{
1106//mDap->SetBlocked();
1107SetState();
1108PIWindow::Show();
1109return;
1110}
1111
1112/* --Methode-- */
1113void PIAFitterWind::SetState(void)
1114// Ecriture et positionnement des ckb et txt
1115{
1116int npar = mFitter->mNPar;
1117
1118// Les check-boxes
1119ckb[0]->SetState(mFitter->mOpt.okfun);
1120ckb[1]->SetState(mFitter->mOpt.okres);
1121if(npar>0)
1122 for(int i=0;i<npar;i++) ckb[2+i]->SetState((mFitter->mFix(i)>=1.)?true:false);
1123
1124// Les champs textes
1125char str[128]; string dum;
1126sprintf(str,"%f",ErrValue); dum=str; txt[0]->SetText(dum);
1127sprintf(str,"%d,%d",mFitter->mOpt.lp,mFitter->mOpt.lpg); dum=str; txt[1]->SetText(dum);
1128sprintf(str,"%f",mFitter->mOpt.stc2); dum=str; txt[2]->SetText(dum);
1129sprintf(str,"%d",mFitter->mOpt.nstep); dum=str; txt[3]->SetText(dum);
1130sprintf(str,"%d",mFitter->mOpt.i1); dum=str; txt[4]->SetText(dum);
1131sprintf(str,"%d",mFitter->mOpt.i2); dum=str; txt[5]->SetText(dum);
1132sprintf(str,"%d",mFitter->mOpt.j1); dum=str; txt[6]->SetText(dum);
1133sprintf(str,"%d",mFitter->mOpt.j2); dum=str; txt[7]->SetText(dum);
1134sprintf(str,"%f",mFitter->mOpt.xc); dum=str; txt[8]->SetText(dum);
1135sprintf(str,"%f",mFitter->mOpt.yc); dum=str; txt[9]->SetText(dum);
1136if(npar>0) {
1137 for(int i=0;i<npar;i++) {
1138 sprintf(str,"%f",mFitter->mPar(i)); dum=str; txt[10+4*i+0]->SetText(dum);
1139 sprintf(str,"%f",mFitter->mStep(i)); dum=str; txt[10+4*i+1]->SetText(dum);
1140 sprintf(str,"%f",mFitter->mMin(i)); dum=str; txt[10+4*i+2]->SetText(dum);
1141 sprintf(str,"%f",mFitter->mMax(i)); dum=str; txt[10+4*i+3]->SetText(dum);
1142 }
1143}
1144
1145return;
1146}
1147
1148/* --Methode-- */
1149void PIAFitterWind::Process(PIMessage msg, PIMsgHandler* sender, void* data)
1150{
1151int lp = (mFitter->mOpt.lp>2)?1:0;;
1152int npar = mFitter->mNPar;
1153
1154// Decodage des messages, init et dismiss
1155msg = UserMsg(msg);
1156if(msg ==777) {
1157 if(lp) cout<<"Dismiss"<<endl;
1158 //mDap->SetReady();
1159 Hide();
1160 return;
1161}
1162
1163// L'objet existe t-il encore?
1164NamedObjMgr omg;
1165if(omg.GetObj(mFitter->mNObj) == NULL) {
1166 cout<<"PIAFitterWind::Process Error , Pas d'objet de nom "
1167 <<mFitter->mNObj<<endl;
1168 mFitter->ResetDPointer(); return;
1169}
1170
1171// On recupere les champs textes
1172string dum; char str[128];
1173// Attention gestion du type d'erreur
1174dum=txt[0]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1175 sscanf(str,"%lf",&ErrValue); if(ErrValue<=0.) ErrValue=1.;
1176 if(mFitter->mOpt.err_e>0. ) mFitter->mOpt.err_e = ErrValue;
1177 if(mFitter->mOpt.err_E>0. ) mFitter->mOpt.err_E = ErrValue;
1178dum=txt[1]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1179 sscanf(str,"%d,%d",&(mFitter->mOpt.lp),&(mFitter->mOpt.lpg));
1180dum=txt[2]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1181 sscanf(str,"%lf",&(mFitter->mOpt.stc2));
1182dum=txt[3]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1183 sscanf(str,"%d",&(mFitter->mOpt.nstep));
1184// Attention si les valeurs de mOpt.i1/.i2/.j1/.j2 ont ete changees
1185dum=txt[4]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1186 sscanf(str,"%d",&(mFitter->mOpt.i1));
1187dum=txt[5]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1188 sscanf(str,"%d",&(mFitter->mOpt.i2));
1189dum=txt[6]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1190 sscanf(str,"%d",&(mFitter->mOpt.j1));
1191dum=txt[7]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1192 sscanf(str,"%d",&(mFitter->mOpt.j2));
1193 if(saveI1!=mFitter->mOpt.i1 || saveI2!=mFitter->mOpt.i2 ||
1194 saveJ1!=mFitter->mOpt.j1 || saveJ2!=mFitter->mOpt.j2 ) {
1195 if(lp)
1196 cout<<"Range changed: I("<<saveI1<<","<<saveI2
1197 <<") -> I("<<mFitter->mOpt.i1<<","<<mFitter->mOpt.i2
1198 <<") ou J("<<saveJ1<<","<<saveJ2
1199 <<") -> J("<<mFitter->mOpt.j1<<","<<mFitter->mOpt.j2
1200 <<")"<<endl;
1201 ReFillGData = true;
1202 saveI1=mFitter->mOpt.i1; saveI2=mFitter->mOpt.i2;
1203 saveJ1=mFitter->mOpt.j1; saveJ2=mFitter->mOpt.j2;
1204 }
1205dum=txt[8]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1206 sscanf(str,"%lf",&(mFitter->mOpt.xc));
1207dum=txt[9]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1208 sscanf(str,"%lf",&(mFitter->mOpt.yc));
1209
1210if(npar>0) {
1211 for(int i=0;i<npar;i++) {
1212 dum=txt[10+4*i+0]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1213 sscanf(str,"%lf",&(mFitter->mPar(i)));
1214 dum=txt[10+4*i+1]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1215 sscanf(str,"%lf",&(mFitter->mStep(i)));
1216 dum=txt[10+4*i+2]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1217 sscanf(str,"%lf",&(mFitter->mMin(i)));
1218 dum=txt[10+4*i+3]->GetText(); strcpy(str,dum.c_str()); strip(str,'B',' ');
1219 sscanf(str,"%lf",&(mFitter->mMax(i)));
1220 }
1221}
1222
1223// Decodage des messages, suite
1224if(msg ==666) {
1225 if(lp) cout<<"Check"<<endl;
1226 mFitter->CheckOptions();
1227 if(ReFillGData) mFitter->FillGData();
1228 ReFillGData = false;
1229 if(mFitter->mFit) mFitter->mFit->PrintStatus();
1230}
1231else if(msg ==111) {
1232 if(lp) cout<<"Update from last"<<endl;
1233 mFitter->ResetOptions();
1234 mFitter->mOpt = mFitter->mOptSave;
1235 mFitter->ReSetParam();
1236 mFitter->InitParFromLastFit();
1237 ReFillGData = true;
1238}
1239else if(msg ==222) {
1240 if(lp) cout<<"Default options/par"<<endl;
1241 mFitter->ResetOptions();
1242 mFitter->ReSetParam();
1243 ReFillGData = true;
1244}
1245else if(msg ==333) {
1246 if(lp) cout<<"Do the fit"<<endl;
1247 mFitter->CheckOptions();
1248 if(ReFillGData) mFitter->FillGData();
1249 ReFillGData = false;
1250 int rc = mFitter->DoFit();
1251 if(rc>=0) mFitter->FitFunRes();
1252}
1253else if(msg==1001) {
1254 if(ckb[0]->GetState()) {
1255 if(lp) cout<<"Gen fitted function"<<endl;
1256 mFitter->mOpt.okfun = true;
1257 } else if(!ckb[0]->GetState()) {
1258 if(lp) cout<<"No Gen fitted function"<<endl;
1259 mFitter->mOpt.okfun = false;
1260 }
1261}
1262else if(msg==1002) {
1263 if(ckb[1]->GetState()) {
1264 if(lp) cout<<"Gen fitted residus"<<endl;
1265 mFitter->mOpt.okres = true;
1266 } else if(!ckb[1]->GetState()) {
1267 if(lp) cout<<"No Gen fitted residus"<<endl;
1268 mFitter->mOpt.okres = false;
1269 }
1270}
1271else if(2101<=msg && msg<=2103) {
1272 mFitter->mOpt.err_e = mFitter->mOpt.err_E = -1.;
1273 if(msg==2102) mFitter->mOpt.err_e = ErrValue;
1274 else if(msg==2103) mFitter->mOpt.err_E = ErrValue;
1275 if(lp) cout<<"Err_cste="<<mFitter->mOpt.err_e
1276 <<" Err_sqrt="<<mFitter->mOpt.err_E
1277 <<" ErrValue="<<ErrValue<<endl;
1278 ReFillGData = true;
1279}
1280else if(2001<=msg && msg<=2003) {
1281 mFitter->mOpt.polcx = msg-2001;
1282 if(lp) cout<<"Centrage X polcx="<<mFitter->mOpt.polcx
1283 <<" xc="<<mFitter->mOpt.xc<<endl;
1284 ReFillGData = true;
1285}
1286else if(2011<=msg && msg<=2013) {
1287 mFitter->mOpt.polcx = msg-2011;
1288 if(lp) cout<<"Centrage Y polcy="<<mFitter->mOpt.polcy
1289 <<" yc="<<mFitter->mOpt.yc<<endl;
1290 ReFillGData = true;
1291}
1292else if(msg>=1500 && msg<2000) {
1293 int ip = msg-1500;
1294 if(ckb[2+ip]->GetState()) {
1295 if(lp) cout<<"Fix Param "<<ip<<endl;
1296 mFitter->mFix(ip) = 1.001;
1297 } else if(!ckb[2+ip]->GetState()) {
1298 if(lp) cout<<"Not Fix Param "<<ip<<endl;
1299 mFitter->mFix(ip) = 0.;
1300 }
1301}
1302
1303mFitter->CheckOptions();
1304SetState();
1305if(mFitter->mOpt.lp>1) mFitter->PrintOptions();
1306
1307return;
1308}
Note: See TracBrowser for help on using the repository browser.