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

Last change on this file since 3455 was 2815, checked in by ansari, 20 years ago

Suite correction PIAFitter en multithread - ajout mutex d'execution pour PIAFitter , cmv+Reza 8 Juillet 2005

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