- Timestamp:
- Aug 6, 1999, 7:16:22 PM (26 years ago)
- Location:
- trunk/SophyaPI/PIext
- Files:
-
- 2 added
- 8 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/SophyaPI/PIext/basexecut.cc
r357 r361 566 566 } 567 567 568 569 // Fit 1D sur objets 1D. Egalement Fit 2D sur objets 2D.570 else if (kw == "fit") {571 if (tokens.size() < 2) {572 cout <<"Usage:fit nomobj func \n"573 <<" [p:p1,...,pn s:s1,...,sn m:m1,...,mn M:M1,...,Mn o:... o:...]\n";574 return(0);575 }576 string p=""; string s=""; string m=""; string M=""; string O="";577 if (tokens.size()>2)578 for(int ip=2;ip<tokens.size();ip++) {579 if(tokens[ip].length()<=2) continue;580 const char *c = tokens[ip].c_str();581 if(c[1]!=':') continue;582 if(c[0]=='p') p=c+2;583 else if(c[0]=='s') s=c+2;584 else if(c[0]=='m') m=c+2;585 else if(c[0]=='M') M=c+2;586 else if(c[0]=='o') {O += ","; O += c+2;}587 }588 srvo->Fit12D(tokens[0],tokens[1],p,s,m,M,O);589 }590 568 591 569 … … 918 896 mpiac->RegisterCommand(kw, usage, this, "Expr. Plotting"); 919 897 920 kw = "fit";921 usage = "Fitting function to DataObjects (Histo, Histo2D, Vector, ...)";922 usage += "\n Usage: fit nomobj func [Options]";923 usage += "\n [p:p1,...,pn s:s1,...,sn m:m1,...,mn M:M1,...,Mn o:... o:...]";924 mpiac->RegisterCommand(kw, usage, this, "Expr. Plotting");925 926 898 } 927 899 -
trunk/SophyaPI/PIext/nomhistadapter.cc
r344 r361 298 298 } 299 299 300 300 //------------------------------------------------------------------------- 301 // Class Adaptateur d'objet (Pour NamedObjMgr) d'objet XNTuple 302 //------------------------------------------------------------------------- 303 304 /* --Methode-- */ 305 NOMAdapter_XNTuple::NOMAdapter_XNTuple(XNTuple* o) 306 : NObjMgrAdapter(o) 307 { 308 mNt = o; 309 } 310 311 /* --Methode-- */ 312 NOMAdapter_XNTuple::~NOMAdapter_XNTuple() 313 { 314 } 315 316 /* --Methode-- */ 317 NObjMgrAdapter* NOMAdapter_XNTuple::Clone(AnyDataObj* o) 318 { 319 XNTuple* nt = dynamic_cast<XNTuple *>(o); 320 if (nt) return ( new NOMAdapter_XNTuple(nt) ); 321 return ( new NObjMgrAdapter(o) ); 322 } 323 324 325 /* --Methode-- */ 326 void NOMAdapter_XNTuple::SavePPF(POutPersist& pos, string const & nom) 327 { 328 #ifdef SANS_EVOLPLANCK 329 // PEIDA-EROS L'histo est lui-meme PPersist 330 string tag = nom; // A cause de const 331 mNt->Write(pos,0,tag); 332 #else 333 string s = typeid(*mObj).name(); 334 cout << "NOMAdapter_XNTuple::SavePPF() - Error : Not supported for " << s << endl; 335 #endif 336 } 337 338 /* --Methode-- */ 339 void NOMAdapter_XNTuple::Print(ostream& os) 340 { 341 // os << mNt->Info(); 342 mNt->Show(os); 343 } 344 345 346 /* --Methode-- */ 347 NTupleInterface* NOMAdapter_XNTuple::GetNTupleInterface(bool& adel) 348 { 349 adel = false; 350 return(mNt); 351 } 352 353 -
trunk/SophyaPI/PIext/nomhistadapter.h
r344 r361 12 12 #include "hisprof.h" 13 13 #include "ntuple.h" 14 #include "xntuple.h" 14 15 15 16 //------------------------------------------------------------------------- … … 95 96 class NOMAdapter_NTuple : public NObjMgrAdapter { 96 97 public: 97 NOMAdapter_NTuple(NTuple* h= NULL);98 NOMAdapter_NTuple(NTuple* nt = NULL); 98 99 virtual ~NOMAdapter_NTuple(); 99 100 … … 111 112 }; 112 113 114 //------------------------------------------------------------------------- 115 // Class Adaptateur d'objet (Pour NamedObjMgr) d'objet XNTuple 116 //------------------------------------------------------------------------- 117 118 class NOMAdapter_XNTuple : public NObjMgrAdapter { 119 public: 120 NOMAdapter_XNTuple(XNTuple* nt = NULL); 121 virtual ~NOMAdapter_XNTuple(); 122 123 virtual NObjMgrAdapter* Clone(AnyDataObj* o); 124 125 // virtual void ReadFits(string const & flnm); 126 // virtual void SaveFits(string const & flnm); 127 virtual void SavePPF(POutPersist& s, string const & nom); 128 129 virtual void Print(ostream& os); 130 virtual NTupleInterface* GetNTupleInterface(bool& adel); 131 132 protected: 133 XNTuple* mNt; 134 }; 135 113 136 114 137 -
trunk/SophyaPI/PIext/piacmd.cc
r357 r361 13 13 #include "pistdimgapp.h" 14 14 #include "nobjmgr.h" 15 #include "piafitting.h" 15 16 16 17 #include PISTDWDG_H … … 309 310 310 311 basexec = new PIABaseExecutor(this, omg, app); 312 fitexec = new PIAFitter(this, app); 311 313 AddInterpreter(this); 312 314 curcmdi = this; … … 327 329 delete helpwin; 328 330 if (curpiacmd == this) curpiacmd = NULL; 331 delete basexec; 332 delete fitexec; 329 333 } 330 334 -
trunk/SophyaPI/PIext/piacmd.h
r351 r361 97 97 CmdInterpreter* curcmdi; 98 98 CmdExecutor* basexec; 99 CmdExecutor* fitexec; 99 100 100 101 // Pour enregistrer la liste de commandes et leurs executeurs et le help -
trunk/SophyaPI/PIext/piinit.cc
r339 r361 23 23 serv->RegisterClass(new Histo2D, new NOMAdapter_Histo2D ); 24 24 serv->RegisterClass(new NTuple, new NOMAdapter_NTuple ); 25 serv->RegisterClass(new XNTuple, new NOMAdapter_XNTuple ); 25 26 26 27 serv->RegisterClass(new GeneralFitData, new NOMAdapter_GeneralFitData ); -
trunk/SophyaPI/PIext/servnobjm.cc
r357 r361 961 961 delete nt; 962 962 mOmg->AddObj(gfd, nomgfd); 963 return;964 }965 966 967 ///////////////////// Fit 1D et 2D //////////////////////////968 /* --Function static propre aux routines de fit 1D et 2D-- cmv 13/10/98 */969 struct DFOptions {970 int okres, okfun;971 int polcx,polcy; double xc,yc;972 double err_e, err_E;973 double stc2;974 int nstep;975 int lp,lpg;976 int i1,i2,j1,j2;977 };978 typedef struct DFOptions DFOPTIONS;979 static void DecodeFitsOptions(string par,string step,string min,string max,string opt980 ,Vector& Par,Vector& Step,Vector& Min,Vector& Max,DFOPTIONS& O);981 void DecodeFitsOptions(string par,string step,string min,string max,string opt982 ,Vector& Par,Vector& Step,Vector& Min,Vector& Max,DFOPTIONS& O)983 //| Pour decoder les "string" et remplir les vecteurs du fit (cf commentaires dans Fit1D)984 {985 // set des vecteurs et decodage des string correspondantes986 int NParMax = 100;987 Par.Realloc(NParMax); Step.Realloc(NParMax);988 Min.Realloc(NParMax); Max.Realloc(NParMax);989 {990 Vector* v=NULL; string* s=NULL;991 {for(int i=0;i<NParMax;i++) {Par(i)=0.; Step(i)=1.; Min(i)=1.; Max(i)=-1.;}}992 for(int j=0;j<4;j++) {993 if(j==0) {v=&Par; s=∥}994 else if(j==1) {v=&Step; s=&step;}995 else if(j==2) {v=&Min; s=&min;}996 else if(j==3) {v=&Max; s=&max;}997 if(s->length()>0) *s += ",";998 for(int i=0;i<NParMax;i++) {999 if(s->length()<=0) break;1000 sscanf(s->c_str(),"%lf",&(*v)(i));1001 size_t p = s->find_first_of(',') + 1;1002 if(p>=s->length()) *s = ""; else *s = s->substr(p);1003 }1004 }1005 }1006 1007 // Decodage de options de opt1008 O.okres = O.okfun = 0;1009 O.polcx = O.polcy = 0;1010 O.xc = O.yc = 0.;1011 O.stc2 = 1.e-3;1012 O.nstep = 100;1013 O.err_e = O.err_E = -1.;1014 O.lp = 1; O.lpg = 0;1015 O.i1 = O.j1 = O.i2 = O.j2 = -1;1016 1017 if(opt.length()<=0) return;1018 opt = "," + opt + ",";1019 1020 if(strstr(opt.c_str(),",r,")) O.okres = 1; // residus1021 if(strstr(opt.c_str(),",f,")) O.okfun = 1; // fonction fittee1022 if(strstr(opt.c_str(),",x")) { // Demande de centrage (fit X=x-xc)1023 O.polcx = 2; // Le centrage est calcule automatiquement1024 size_t p = opt.find(",x");1025 size_t q = opt.find_first_of(',',p+1);1026 string dum = opt.substr(p,q-p);1027 if(dum.length()>2) {1028 sscanf(dum.c_str(),",x%lf",&O.xc);1029 O.polcx = 1; // Le centrage est fixe par la valeur lue1030 }1031 }1032 if(strstr(opt.c_str(),",y")) { // Demande de centrage (fit Y=y-yc)1033 O.polcy = 2; // Le centrage est calcule automatiquement1034 size_t p = opt.find(",y");1035 size_t q = opt.find_first_of(',',p+1);1036 string dum = opt.substr(p,q-p);1037 if(dum.length()>2) {1038 sscanf(dum.c_str(),",y%lf",&O.yc);1039 O.polcy = 1; // Le centrage est fixe par la valeur lue1040 }1041 }1042 if(strstr(opt.c_str(),",E")) { // Erreurs imposees a "sqrt(val)" ou "aa.b*sqrt(val)"1043 size_t p = opt.find(",E");1044 size_t q = opt.find_first_of(',',p+1);1045 string dum = opt.substr(p,q-p);1046 if(dum.length()>2) sscanf(dum.c_str(),",E%lf",&O.err_E);1047 if(O.err_E<=0.) O.err_E = 1.;1048 O.err_e=-1.;1049 }1050 if(strstr(opt.c_str(),",e")) { // Erreurs imposees a "1" ou "aa.b"1051 size_t p = opt.find(",e");1052 size_t q = opt.find_first_of(',',p+1);1053 string dum = opt.substr(p,q-p);1054 if(dum.length()>2) sscanf(dum.c_str(),",e%lf",&O.err_e);1055 if(O.err_e<=0.) O.err_e = 1.;1056 O.err_E=-1.;1057 }1058 if(strstr(opt.c_str(),",X")) { // Valeur du StopChi21059 size_t p = opt.find(",X");1060 size_t q = opt.find_first_of(',',p+1);1061 string dum = opt.substr(p,q-p);1062 if(dum.length()>2) sscanf(dum.c_str(),",X%lf",&O.stc2);1063 if(O.stc2<=0.) O.stc2 = 1.e-3;1064 }1065 if(strstr(opt.c_str(),",N")) { // Nombre maximum d'iterations1066 size_t p = opt.find(",N");1067 size_t q = opt.find_first_of(',',p+1);1068 string dum = opt.substr(p,q-p);1069 if(dum.length()>2) sscanf(dum.c_str(),",N%d",&O.nstep);1070 if(O.nstep<2) O.nstep = 100;1071 }1072 if(strstr(opt.c_str(),",l")) { // niveau de print1073 size_t p = opt.find(",l");1074 size_t q = opt.find_first_of(',',p+1);1075 string dum = opt.substr(p,q-p);1076 float ab;1077 if(dum.length()>2) sscanf(dum.c_str(),",l%f",&ab);1078 if(ab<0) ab = 0.;1079 O.lp = (int) ab; O.lpg = int(10.*(ab-(float)O.lp+0.01));1080 }1081 if(strstr(opt.c_str(),",I")) { // intervalle de fit selon X1082 size_t p = opt.find(",I");1083 size_t q = opt.find_first_of(',',p+1);1084 string dum = opt.substr(p,q-p);1085 if(dum.length()>2) sscanf(dum.c_str(),",I%d/%d",&O.i1,&O.i2);1086 }1087 if(strstr(opt.c_str(),",J")) { // intervalle de fit selon Y1088 size_t p = opt.find(",J");1089 size_t q = opt.find_first_of(',',p+1);1090 string dum = opt.substr(p,q-p);1091 if(dum.length()>2) sscanf(dum.c_str(),",J%d/%d",&O.j1,&O.j2);1092 }1093 return;1094 }1095 1096 /* --Methode-- cmv 13/10/98 */1097 void Services2NObjMgr::Fit12D(string& nom, string& func,1098 string par,string step,string min,string max,1099 string opt)1100 //| --------------- Fit d'objets a 1 et 2 dimensions ---------------1101 //| nom : nom de l'objet qui peut etre:1102 //| fit-1D: Vector,Histo1D,HProf ou GeneraFitData(1D)1103 //| fit-2D: Matrix,Histo2D,Image<T> ou GeneraFitData(2D)1104 //| func : pnn : fit polynome degre nn avec classe Poly (lineaire) 1D ou 2D1105 //| : Pnn : fit polynome degre nn avec GeneralFit (non-lineaire) 1D ou 2D1106 //| : gnn : fit gaussienne (hauteur) + polynome de degre nn 1D1107 //| : g : fit gaussienne (hauteur) 1D1108 //| : enn : fit exponentielle + polynome de degre nn 1D1109 //| : e : fit exponentielle 1D1110 //| : Gnn : fit gaussienne (volume) + polynome de degre nn 1D1111 //| : G : fit gaussienne (volume) 1D1112 //| : : fit gaussienne+fond (volume) 2D1113 //| : Gi : fit gaussienne+fond integree (volume) 2D1114 //| : d : fit DL de gaussienne+fond (volume) 2D1115 //| : di : fit DL de gaussienne+fond integree (volume) 2D1116 //| : D : fit DL de gaussienne+fond avec coeff variable p6 (volume) 2D1117 //| : Di : fit DL de gaussienne+fond integree avec coeff variable p6 (volume) 2D1118 //| : M : fit Moffat+fond (expos=p6) (volume) 2D1119 //| : Mi : fit Moffat+fond integree (expos=p6) (volume) 2D1120 //| par : p1,...,pn : valeur d'initialisation des parametres (def=0)1121 //| step : s1,...,sn : valeur des steps de depart (def=1)1122 //| min : m1,...,mn : valeur des minima (def=1)1123 //| max : M1,...,Mn : valeur des maxima (def=-1) (max<=min : pas de limite)1124 //| opt : options "Eaa.b,eaa.b,f,r,caa.b,Xaa.b"1125 //| f : generation d'un Objet identique contenant la fonction fittee1126 //| r : generation d'un Objet identique contenant les residus1127 //| Xaa.b : aa.b valeur du DXi2 d'arret (def=1.e-3)1128 //| Naa : aa nombre maximum d'iterations (def=100)1129 //| la.b : niveau "a.b" de print: a=niveau de print Fit1/2D1130 //| b=niveau de debug GeneralFit1131 //| Ii1/i2 numeros des bins X de l'histos utilises pour le fit [i1,i2]1132 //|2D Jj1/j2 numeros des bins Y de l'histos utilises pour le fit [j1,j2]1133 //| - L'erreur est celle associee a l'objet (si elle existe),1134 //| elle est mise a 1 sinon, sauf si E... ou e... est precise:1135 //| Eaa.b : si |val|>=1 erreur = aa.b*sqrt(|val|)1136 //| si |val|<1 erreur = aa.b1137 //| si aa.b <=0 alors aa.b=1.01138 //| E seul est equivalent a E1.01139 //| eaa.b : erreur = aa.b1140 //| si aa.b <=0 alors aa.b=1.01141 //| e seul est equivalent a e1.01142 //| xaa.b : demande de centrage: on fit x-aa.b au lieu de x)1143 //| x : demande de centrage: on fit x-xc au lieu de x1144 //| avec xc=abscisse du milieu de l'histogramme1145 //| Actif pour exp+poly 1D, poly 1D1146 //| pour gauss+poly 1D, xc est le centre de la gaussienne.1147 //|2D yaa.b et y : idem "xaa.b et x" mais pour y1148 {1149 AnyDataObj* obj=mOmg->GetObj(nom);1150 if (obj == NULL) {1151 cout<<"Services2NObjMgr::Fit12D() Error , Pas d'objet de nom "<<nom<<endl;1152 return;1153 }1154 if (!mImgapp) return;1155 if(func.length()<=0)1156 {cout<<"Services2NObjMgr::Fit12D() Donnez un nom de fonction a fitter."<<endl;1157 return;}1158 string ctyp = typeid(*obj).name();1159 1160 int ndim = 0, nbinx=0, nbiny=0, ndata = 0;1161 Vector* v = NULL; Histo* h = NULL;1162 Matrix* m = NULL; Histo2D* h2 = NULL; RzImage* im = NULL;1163 GeneralFitData* g = NULL;1164 1165 // 1D1166 if (typeid(*obj) == typeid(Vector)) {1167 ndim = 1;1168 v = (Vector*) obj; nbinx = v->NElts(); nbiny = 1;1169 }1170 else if ( (typeid(*obj) == typeid(HProf)) || (typeid(*obj) == typeid(Histo)) ) {1171 ndim = 1;1172 h = (Histo*) obj; nbinx = h->NBins(); nbiny = 1;1173 }1174 else if (typeid(*obj) == typeid(Matrix)) {1175 ndim = 2;1176 m = (Matrix*) obj; nbinx = m->NCol(); nbiny = m->NRows();1177 }1178 else if (typeid(*obj) == typeid(Histo2D)) {1179 ndim = 2;1180 h2 = (Histo2D*) obj; nbinx = h2->NBinX(); nbiny = h2->NBinY();1181 }1182 else if (typeid(*obj) == typeid(GeneralFitData)) {1183 g = (GeneralFitData*) obj; nbinx = g->NData(); nbiny = 1;1184 if( g->NVar()==1) ndim = 1;1185 else if(g->NVar()==2) ndim = 2;1186 else {1187 cout<<"GeneralFitData ne peut avoir que 1 ou 2 variables d'abscisse: "1188 <<((GeneralFitData*) obj)->NVar()<<endl; return; }1189 }1190 else if (dynamic_cast<RzImage*>(obj)) {1191 ndim = 2;1192 im = (RzImage*) obj; nbinx = im->XSize(); nbiny = im->YSize();1193 }1194 else {1195 cout<<"Services2NObjMgr::Fit12D() Error , Objet n'est pas un "1196 <<"Histo1D/HProf/Vector/Histo2D/Image/Matrix/GeneralFitData "<<ctyp<<endl;1197 return;1198 }1199 1200 ndata = nbinx*nbiny;1201 if(ndata<=0)1202 {cout<<"L'objet a "<<nbinx<<","<<nbiny<<" bins ("<<ndata<<")"<<endl; return;}1203 1204 // Decodage des options et des parametres, mise en forme1205 Vector Par(1); Vector Step(1); Vector Min(1); Vector Max(1); DFOPTIONS O;1206 DecodeFitsOptions(par,step,min,max,opt,Par,Step,Min,Max,O);1207 O.i1 = (O.i1<0||O.i1>=nbinx)? 0: O.i1;1208 O.i2 = (O.i2<0||O.i2>=nbinx||O.i2<O.i1)? nbinx-1: O.i2;1209 if(ndim>=2) {1210 O.j1 = (O.j1<0||O.j1>=nbiny)? 0: O.j1;1211 O.j2 = (O.j2<0||O.j2>=nbiny||O.j2<O.j1)? nbiny-1: O.j2;1212 } else O.j2 = O.j1 = 0;1213 if(O.polcx==2) {1214 if(v||m) O.xc = (O.i2-O.i1+1)/2.;1215 else if(h) O.xc = (h->XMin()+h->XMax())/2.;1216 else if(h2) O.xc = (h2->XMin()+h2->XMax())/2.;1217 else if(g) {double mini,maxi; g->GetMinMax(2,mini,maxi); O.xc=(mini+maxi)/2.;}1218 else if(im) {O.xc = im->XOrg() * im->XPxSize()*(O.i2-O.i1+1)/2.;}1219 }1220 if(O.polcy==2 && ndim>=2) {1221 if(m) O.yc = (O.j2-O.j1+1)/2.;1222 if(h2) O.yc = (h2->YMin()+h2->YMax())/2.;1223 if(g) {double mini,maxi; g->GetMinMax(12,mini,maxi); O.yc=(mini+maxi)/2.;}1224 if(im) {O.yc = im->YOrg() * im->YPxSize()*(O.j2-O.j1+1)/2.;}1225 }1226 if(O.lp>0)1227 cout<<"Fit["<<nbinx<<","<<nbiny<<"] ("<<ndata<<") dim="<<ndim<<":"1228 <<" Int=["<<O.i1<<","<<O.i2<<"],["<<O.j1<<","<<O.j2<<"]"<<endl1229 <<" Cent="<<O.polcx<<","<<O.polcy<<","<<O.xc<<"+x"<<","<<O.yc<<"+y"1230 <<" TypE="<<O.err_e<<","<<O.err_E1231 <<" StpX2="<<O.stc2<<" Nstep="<<O.nstep1232 <<" lp,lpg="<<O.lp<<","<<O.lpg<<endl;1233 1234 ///////////////////////////////////1235 // Remplissage de GeneralFitData //1236 ///////////////////////////////////1237 GeneralFitData mydata(ndim,ndata,0);1238 {for(int i=O.i1;i<=O.i2;i++) for(int j=O.j1;j<=O.j2;j++) {1239 double x,y,f,e;1240 1241 if(v)1242 {x= (double) i; f=(*v)(i); e=1.;}1243 else if(h)1244 {x=h->BinCenter(i); f=(*h)(i); e=(h->HasErrors())?h->Error(i):1.;}1245 else if(m)1246 {x=(double) i; y=(double) j; f=(*m)(j,i); e=1.;}1247 else if(h2)1248 {float xf,yf; h2->BinCenter(i,j,xf,yf); x=(double)xf; y=(double)yf;1249 f=(*h2)(i,j); e=(h2->HasErrors())?h2->Error(i,j):1.;}1250 else if(im)1251 {x=im->XOrg()+(i+0.5)*im->XPxSize(); y=im->YOrg()+(j+0.5)*im->YPxSize();1252 f=im->DValue(i,j); e=1.;}1253 else if(g&&ndim==1) {x= g->X(i); f=g->Val(i); e=g->EVal(i);}1254 else if(g&&ndim==2) {x= g->X(i); y= g->Y(i); f=g->Val(i); e=g->EVal(i);}1255 else x=y=f=e=0.;1256 1257 // Gestion des erreurs a utiliser1258 if(O.err_e>0.) e=O.err_e;1259 else if(O.err_E>0.) {e=(f<-1.||f>1.)?O.err_E*sqrt(fabs(f)):O.err_E;}1260 1261 // Remplissage de generalfit1262 if(func[0]=='p') {x -= O.xc; if(ndim>=2) y -= O.yc;}1263 if(ndim==1) mydata.AddData1(x,f,e);1264 else if(ndim==2) mydata.AddData2(x,y,f,e);1265 }}1266 if(mydata.NData()<=0)1267 {cout<<"Pas de donnees dans GeneralFitData: "<<mydata.NData()<<endl;1268 return;}1269 if(O.lpg>1) {1270 mydata.PrintStatus();1271 mydata.PrintData(0);1272 mydata.PrintData(mydata.NData()-1);1273 }1274 1275 ////////////////////////////////////////////1276 // Identification de la fonction a fitter //1277 ////////////////////////////////////////////1278 GeneralFunction* myfunc = NULL;1279 if(func[0]=='p' && ndim==1) {1280 // Fit de polynome sans passer par les GeneralFit1281 int degre = 0;1282 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1283 cout<<"Fit (lineaire) 1D polynome de degre "<<degre<<endl;1284 Poly p1(0);1285 double c2rl = mydata.PolFit(0,p1,degre);1286 cout<<"C2r_lineaire = "<<c2rl<<endl;1287 if(O.lp>0) cout<<p1<<endl;1288 return;1289 1290 } else if(func[0]=='P' && ndim==1) {1291 // Fit de polynome1292 int degre = 0;1293 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1294 cout<<"Fit polynome 1D de degre "<<degre<<endl;1295 Polyn1D* myf = new Polyn1D(degre,O.xc);1296 myfunc = myf;1297 1298 } else if(func[0]=='e' && ndim==1) {1299 // Fit d'exponentielle1300 int degre =-1;1301 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1302 cout<<"Fit d'exponentielle+polynome 1D de degre "<<degre<<endl;1303 Exp1DPol* myf;1304 if(degre>=0) myf = new Exp1DPol((unsigned int)degre,O.xc);1305 else myf = new Exp1DPol(O.xc);1306 myfunc = myf;1307 1308 } else if(func[0]=='g' && ndim==1) {1309 // Fit de gaussienne en hauteur1310 int degre =-1;1311 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1312 cout<<"Fit de Gaussienne_en_hauteur+polynome 1D de degre "<<degre<<endl;1313 Gauss1DPol* myf;1314 if(degre>=0) myf = new Gauss1DPol((unsigned int)degre,((O.polcx)?true:false));1315 else { bool bfg = (O.polcx)?true:false; myf = new Gauss1DPol(bfg); }1316 myfunc = myf;1317 1318 } else if(func[0]=='G' && ndim==1) {1319 // Fit de gaussienne en volume1320 int degre =-1;1321 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1322 cout<<"Fit de Gaussienne_en_volume+polynome 1D de degre "<<degre<<endl;1323 GaussN1DPol* myf;1324 if(degre>=0) myf = new GaussN1DPol((unsigned int)degre,((O.polcx)?true:false));1325 else { bool bfg = (O.polcx)?true:false; myf = new GaussN1DPol(bfg); }1326 myfunc = myf;1327 1328 } else if(func[0]=='p' && ndim==2) {1329 // Fit de polynome 2D sans passer par les GeneralFit1330 int degre = 0;1331 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1332 cout<<"Fit (lineaire) polynome 2D de degre "<<degre<<endl;1333 Poly2 p2(0);1334 double c2rl = mydata.PolFit(0,1,p2,degre);1335 cout<<"C2r_lineaire = "<<c2rl<<endl;1336 if(O.lp>0) cout<<p2<<endl;1337 return;1338 1339 } else if(func[0]=='P' && ndim==2) {1340 // Fit de polynome 2D1341 int degre = 0;1342 if(func.length()>1) sscanf(func.c_str()+1,"%d",°re);1343 cout<<"Fit polynome 2D de degre "<<degre<<endl;1344 Polyn2D* myf = new Polyn2D(degre,O.xc,O.yc);1345 myfunc = myf;1346 1347 } else if(func[0]=='G' && ndim==2) {1348 // Fit de gaussienne+fond en volume1349 int integ = 0;1350 if(func.length()>1) if(func[1]=='i') integ=1;1351 cout<<"Fit de Gaussienne+Fond 2D integ="<<integ<<endl;1352 if(integ) {GauRhInt2D* myf = new GauRhInt2D; myfunc = myf;}1353 else {GauRho2D* myf = new GauRho2D; myfunc = myf;}1354 1355 } else if(func[0]=='d' && ndim==2) {1356 // Fit de DL gaussienne+fond en volume1357 int integ = 0;1358 if(func.length()>1) if(func[1]=='i') integ=1;1359 cout<<"Fit de DL de Gaussienne+Fond 2D integ="<<integ<<endl;1360 if(integ) {GdlRhInt2D* myf = new GdlRhInt2D; myfunc = myf;}1361 else {GdlRho2D* myf = new GdlRho2D; myfunc = myf;}1362 1363 } else if(func[0]=='D' && ndim==2) {1364 // Fit de DL gaussienne+fond avec coeff variable p6 en volume1365 int integ = 0;1366 if(func.length()>1) if(func[1]=='i') integ=1;1367 cout<<"Fit de DL de Gaussienne+Fond avec coeff variable (p6) 2D integ="<<integ<<endl;1368 if(integ) {Gdl1RhInt2D* myf = new Gdl1RhInt2D; myfunc = myf;}1369 else {Gdl1Rho2D* myf = new Gdl1Rho2D; myfunc = myf;}1370 1371 } else if(func[0]=='M' && ndim==2) {1372 // Fit de Moffat+fond (volume)1373 int integ = 0;1374 if(func.length()>1) if(func[1]=='i') integ=1;1375 cout<<"Fit de Moffat+Fond (expos=p6) 2D integ="<<integ<<endl;1376 if(integ) {MofRhInt2D* myf = new MofRhInt2D; myfunc = myf;}1377 else {MofRho2D* myf = new MofRho2D; myfunc = myf;}1378 1379 } else {1380 cout<<"Fonction "<<func<<" inconnue pour la dim "<<ndim<<endl;1381 return;1382 }1383 1384 /////////////////////////1385 // Fit avec generalfit //1386 /////////////////////////1387 if(myfunc->NPar()>Par.NElts())1388 {cout<<"Trop de parametres: "<<myfunc->NPar()<<">"<<Par.NElts()<<endl;1389 if(myfunc) delete myfunc; return;}1390 GeneralFit myfit(myfunc);1391 myfit.SetDebug(O.lpg);1392 myfit.SetData(&mydata);1393 myfit.SetStopChi2(O.stc2);1394 myfit.SetMaxStep(O.nstep);1395 {for(int i=0;i<myfunc->NPar();i++) {1396 char str[10];1397 sprintf(str,"P%d",i);1398 myfit.SetParam(i,str,Par(i),Step(i),Min(i),Max(i));1399 }}1400 if(O.lp>1) myfit.PrintFit();1401 double c2r = -1.;1402 int rcfit = (double) myfit.Fit();1403 if(O.lp>0) myfit.PrintFit();1404 if(rcfit>0) {1405 c2r = myfit.GetChi2Red();1406 cout<<"C2r_Reduit = "<<c2r<<" nstep="<<myfit.GetNStep()<<" rc="<<rcfit<<endl;1407 Vector ParFit(myfunc->NPar());1408 for(int i=0;i<myfunc->NPar();i++) ParFit(i)=myfit.GetParm(i);1409 } else {1410 cout<<"echec Fit, rc = "<<rcfit<<" nstep="<<myfit.GetNStep()<<endl;1411 myfit.PrintFitErr(rcfit);1412 }1413 1414 // Mise a disposition des resultats1415 if(rcfit>=0 && myfunc && (O.okres>0||O.okfun>0)) {1416 string nomres = nom + "res";1417 string nomfun = nom + "fun";1418 if(v) {1419 if(O.okres) {Vector* ob = v->FitResidus(myfit); if(ob) mOmg->AddObj(ob,nomres);}1420 if(O.okfun) {Vector* ob = v->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}1421 } else if(h) {1422 if(O.okres) {Histo* ob = h->FitResidus(myfit); if(ob) mOmg->AddObj(ob,nomres);}1423 if(O.okfun) {Histo* ob = h->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}1424 } else if(m) {1425 if(O.okres) {Matrix* ob = m->FitResidus(myfit); if(ob) mOmg->AddObj(ob,nomres);}1426 if(O.okfun) {Matrix* ob = m->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}1427 } else if(h2) {1428 if(O.okres) {Histo2D* ob = h2->FitResidus(myfit); if(ob) mOmg->AddObj(ob,nomres);}1429 if(O.okfun) {Histo2D* ob = h2->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}1430 } else if(im) {1431 if(O.okres) {RzImage* ob = im->FitResidus(myfit); if(ob) mOmg->AddObj(ob,nomres);}1432 if(O.okfun) {RzImage* ob = im->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}1433 } else if(g) {1434 if(O.okres) {GeneralFitData* ob = g->FitResidus(myfit); if(ob) mOmg->AddObj(ob,nomres);}1435 if(O.okfun) {GeneralFitData* ob = g->FitFunction(myfit); if(ob) mOmg->AddObj(ob,nomfun);}1436 }1437 }1438 1439 // Nettoyage1440 if(myfunc) delete myfunc;1441 963 return; 1442 964 } … … 1932 1454 case ClassId_NTuple : 1933 1455 return("NTuple"); 1456 case ClassId_XNTuple : 1457 return("XNTuple"); 1934 1458 case ClassId_GeneralFitData : 1935 1459 return("GeneralFitData"); -
trunk/SophyaPI/PIext/servnobjm.h
r357 r361 91 91 string const & funcname); 92 92 93 // Methodes de fit - CMV , deux methodes H1 H2 ou 1 seule ??) ...94 virtual void Fit12D(string & nom, string& func,95 string par,string step,string min,string max,string opt);96 97 93 // Calcul d'expressions d'interface NTuple pour les objets 98 94 void ComputeExpressions(NObjMgrAdapter* obja, string& expx, string& expy,
Note:
See TracChangeset
for help on using the changeset viewer.