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

Last change on this file since 1451 was 1268, checked in by ercodmgr, 25 years ago

3ieme argument string& toks pour Execute et ExecuteCommand
preparation gestion TmpDir dans cxxexecutor
suppression de ExecuteCXX ds cxxexecutor et appels
possibilite d'ajouter des userfct.cc ds cxxexecutor

cmv 1/11/00

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