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

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

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