source: Sophya/trunk/SophyaPI/PIext/servnobjm.cc@ 4078

Last change on this file since 4078 was 3572, checked in by cmv, 17 years ago

char* -> const char* pour regler les problemes de deprecated string const... + comparaison unsigned signed + suppression EVOL_PLANCK rz+cmv 07/02/2009

File size: 42.1 KB
Line 
1#include <stdio.h>
2#include <stdlib.h>
3#include <ctype.h>
4
5#include <typeinfo>
6#include <iostream>
7#include <string>
8#include <list>
9#include <map>
10
11#include "sopnamsp.h"
12#include "strutil.h"
13
14#include "nobjmgr.h"
15#include "servnobjm.h"
16#include "nomgadapter.h"
17#include "pistdimgapp.h"
18
19#include "fct1dfit.h"
20#include "fct2dfit.h"
21
22#include "tmatrix.h"
23#include "tvector.h"
24#include "matharr.h"
25#include "pitvmaad.h"
26
27#include "ntuple.h"
28#include "cimage.h"
29
30#include "histos.h"
31#include "histos2.h"
32#include "hisprof.h"
33#include "ntuple.h"
34#include "datatable.h"
35
36#include "piyfxdrw.h"
37#include "pisurfdr.h"
38
39#include "pintuple.h"
40#include "pintup3d.h"
41
42#include "pipodrw.h"
43
44
45
46/* --Methode-- */
47Services2NObjMgr::Services2NObjMgr(string& tmpdir)
48{
49SetTmpDir(tmpdir);
50mImgapp = NULL;
51mOmg = NULL;
52dynlink = NULL;
53}
54
55/* --Methode-- */
56Services2NObjMgr::~Services2NObjMgr()
57{
58CloseDLL();
59if (mOmg) delete mOmg;
60}
61
62/* --Methode-- */
63void Services2NObjMgr::RegisterClass(AnyDataObj* o, NObjMgrAdapter* oa)
64{
65ObjAdaptList::iterator it;
66for(it = objadaplist.begin(); it != objadaplist.end(); it++)
67 if (typeid(*o) == typeid(*((*it).obj)))
68 throw(DuplicateIdExc("Services2NObjMgr::RegisterClass() - Duplicate class"));
69
70dataobj_adapter oba;
71oba.obj = o;
72oba.obja = oa;
73objadaplist.push_back(oba);
74}
75
76/* --Methode-- */
77NObjMgrAdapter* Services2NObjMgr::GetAdapter(AnyDataObj* o)
78{
79ObjAdaptList::iterator it;
80for(it = objadaplist.begin(); it != objadaplist.end(); it++)
81 if (typeid(*o) == typeid(*((*it).obj))) return((*it).obja->Clone(o));
82return(new NObjMgrAdapter(o));
83}
84
85/* --Methode-- */
86void Services2NObjMgr::SetTmpDir(string const & tmpdir)
87{
88TmpDir = tmpdir;
89PDynLinkMgr::SetTmpDir(tmpdir);
90return;
91}
92
93/* --Methode-- */
94void Services2NObjMgr::PlotFunc(string const & expfunc, string & nom, double xmin, double xmax, int np, string dopt)
95{
96FILE *fip;
97string fname = TmpDir + "func1_pia_dl.c";
98string cmd;
99int rc;
100
101if (!mImgapp) return;
102
103// Pour synchronisation d'execution simultanee
104ZSync zs(mutx_dynlink); zs.NOp();
105
106cmd = "rm -f " + fname;
107rc = system(cmd.c_str());
108// printf("PlotFunc_Do> %s (Rc=%d)\n", cmd.c_str(), rc);
109
110if ((fip = fopen(fname.c_str(), "w")) == NULL) {
111 string sn = fname;
112 cout << "Services2NObjMgr/PlotFunc_Error: fopen( " << sn << endl;
113 return;
114 }
115
116// constitution du fichier a compiler
117fputs("#include <math.h> \n", fip);
118fputs("double func1_pia_dl_func(double x) \n{\n", fip);
119fprintf(fip,"return(%s); \n}\n", expfunc.c_str());
120fclose(fip);
121
122string func = "func1_pia_dl_func";
123DlFunctionOfX f = (DlFunctionOfX) LinkFunctionFromFile(fname, func);
124if (!f) return;
125PlotFunc(f, nom, xmin, xmax, np, dopt);
126CloseDLL();
127return;
128}
129
130/* --Methode-- */
131void Services2NObjMgr::PlotFunc2D(string const & expfunc, string & nom, double xmin, double xmax,
132 double ymin, double ymax, int npx, int npy, string dopt)
133{
134FILE *fip;
135string fname = TmpDir + "func2_pia_dl.c";
136string cmd;
137int rc;
138
139if (!mImgapp) return;
140// Pour synchronisation d'execution simultanee
141ZSync zs(mutx_dynlink); zs.NOp();
142
143cmd = "rm " + fname;
144rc = system(cmd.c_str());
145// printf("PlotFunc2D_Do> %s (Rc=%d)\n", cmd.c_str(), rc);
146
147if ((fip = fopen(fname.c_str(), "w")) == NULL) {
148 string sn = fname;
149 cout << "Services2NObjMgr/PlotFunc2D_Error: fopen( " << sn << endl;
150 return;
151 }
152
153// constitution du fichier a compiler
154fputs("#include <math.h> \n", fip);
155fputs("double func2_pia_dl_func(double x, double y) \n{\n", fip);
156fprintf(fip,"return(%s); \n}\n", expfunc.c_str());
157fclose(fip);
158
159string func = "func2_pia_dl_func";
160DlFunctionOfXY f = (DlFunctionOfXY) LinkFunctionFromFile(fname, func);
161if (!f) return;
162PlotFunc2D(f, nom, xmin, xmax, ymin, ymax, npx, npy, dopt);
163CloseDLL();
164return;
165}
166
167/* --Methode-- */
168void Services2NObjMgr::PlotFuncFrCFile(string const & fname, string const & func, string & nom,
169 double xmin, double xmax, int np, string dopt)
170{
171// Pour synchronisation d'execution simultanee
172ZSync zs(mutx_dynlink); zs.NOp();
173DlFunctionOfX f = (DlFunctionOfX) LinkFunctionFromFile(fname, func);
174if (!f) return;
175PlotFunc(f, nom, xmin, xmax, np, dopt);
176CloseDLL();
177return;
178}
179
180/* --Methode-- */
181void Services2NObjMgr::PlotFunc2DFrCFile(string const & fname, string const & func, string & nom,
182 double xmin, double xmax, double ymin, double ymax, int npx, int npy, string dopt)
183{
184// Pour synchronisation d'execution simultanee
185ZSync zs(mutx_dynlink); zs.NOp();
186DlFunctionOfXY f = (DlFunctionOfXY) LinkFunctionFromFile(fname, func);
187if (!f) return;
188PlotFunc2D(f, nom, xmin, xmax, ymin, ymax, npx, npy, dopt);
189CloseDLL();
190return;
191}
192
193/* --Methode-- */
194void Services2NObjMgr::PlotFunc(DlFunctionOfX f, string & nom, double xmin, double xmax, int np, string dopt)
195{
196if (!mImgapp) return;
197
198int k;
199if (np < 1) np = 1;
200if (xmax <= xmin) xmax = xmin+1.;
201Vector* vpy = new Vector(np);
202
203double xx;
204double dx = (xmax-xmin)/np;
205
206try {
207 for(k=0; k<np; k++) { xx = xmin+dx*k; (*vpy)(k) = f(xx); }
208}
209catch ( PException exc ) {
210 fflush(stdout);
211 cout << endl;
212 cerr << endl;
213 cerr << "Services2NObjMgr::PlotFunc() Exception :" << exc.Msg() << endl;
214 delete vpy;
215 vpy = NULL;
216}
217
218if (vpy) {
219 string titre;
220 if (nom.length() < 1) {
221 titre = "Function f(x)";
222 nom = "/autoc/f_x";
223 }
224 else titre = nom;
225
226 MyObjMgr()->AddObj(vpy, nom);
227 if ( (dopt == "nodisp") || (dopt == "nodisplay") ) return;
228
229 P1DArrayAdapter* vya = new POVectorAdapter(vpy, false);
230 vya->DefineXCoordinate(xmin, (xmax-xmin)/np);
231 PIYfXDrawer* dr = new PIYfXDrawer(vya, NULL, true) ;
232 dopt = "thinline " + dopt;
233 int rsid = mImgapp->DispScDrawer(dr, titre, dopt);
234 MyObjMgr()->AddWRsId(nom, rsid);
235}
236return;
237}
238
239/* --Methode-- */
240void Services2NObjMgr::PlotFunc2D(DlFunctionOfXY f, string & nom, double xmin, double xmax, double ymin, double ymax,
241 int npx, int npy, string dopt)
242{
243if (!mImgapp) return;
244
245if (npx < 2) npx = 2;
246if (npy < 2) npy = 2;
247//if (npx > 250) npx = 250;
248//if (npy > 250) npy = 250;
249if (xmax <= xmin) xmax = xmin+1.;
250if (ymax <= ymin) ymax = ymin+1.;
251
252Matrix* mtx = new Matrix(npy, npx);
253
254int i,j;
255double xx, yy;
256double dx = (xmax-xmin)/npx;
257double dy = (ymax-ymin)/npy;
258// printf(" -- DBG -- %d %d , %g %g , %g %g \n", npx, npy, xmin, xmax, ymin, ymax);
259try {
260 for(j=0; j<npy; j++) {
261 yy = ymin+dy*j;
262 for(i=0; i<npx; i++) {
263 xx = xmin+dx*i;
264 (*mtx)(j, i) = f(xx, yy);
265 }
266 }
267}
268catch ( PException exc ) {
269 fflush(stdout);
270 cout << endl;
271 cerr << endl;
272 cerr << "Services2NObjMgr::PlotFunc2D() Exception :" << exc.Msg() << endl;
273 delete mtx; mtx = NULL;
274}
275
276if (mtx) {
277 string titre;
278 if (nom.length() < 1) {
279 titre = "Function f(x,y)";
280 nom = "/autoc/f2d_xy";
281 }
282 else titre = nom;
283
284 MyObjMgr()->AddObj(mtx, nom);
285 if ( (dopt == "nodisp") || (dopt == "nodisplay") ) return;
286
287 int rsid = 0;
288 P2DArrayAdapter* arr = new POMatrixAdapter(mtx, false);
289 arr->DefineXYCoordinates(xmin, ymin, dx, dy);
290 if ( (npx <= 200) && (npy <= 200) ) {
291 PISurfaceDrawer* sdr = new PISurfaceDrawer(arr, true, true, true);
292 rsid = mImgapp->Disp3DDrawer(sdr, titre, dopt);
293 }
294 else rsid = mImgapp->DispImage(arr, titre, dopt);
295 MyObjMgr()->AddWRsId(nom, rsid);
296}
297
298return;
299}
300
301/* --Methode-- */
302void Services2NObjMgr::ExpVal(string expval,string resultvarname)
303{
304 // Fill C-file to be executed
305 FILE *fip;
306 string func = "eval_pia_dl_func";
307 string fname = TmpDir + func; fname += ".c";
308 string cmd = "rm -f " + fname; system(cmd.c_str());
309 if((fip=fopen(fname.c_str(), "w"))==NULL) {
310 cout << "Services2NObjMgr/EvalExp_Error: fopen("<<fname<<")"<<endl;
311 return;
312 }
313 fprintf(fip,"#include <math.h>\n");
314 fprintf(fip,"double %s(double ___dummy_variable___) \n{\n",func.c_str());
315 // Add all variables already declared
316 DVList& varlist = MyObjMgr()->GetVarList();
317 DVList::ValList::const_iterator it;
318 for(it = varlist.Begin(); it != varlist.End(); it++) {
319 double value = (double)((*it).second.elval);
320 string name_var = (*it).first;
321 fprintf(fip,"double %s = %.17f;\n",name_var.c_str(),value);
322 }
323 fprintf(fip,"return %s;\n}\n",expval.c_str());
324 fclose(fip);
325
326 // Dynamically link function
327 DlFunctionOfX f = (DlFunctionOfX) LinkFunctionFromFile(fname,func);
328 if(!f) {
329 cout<<"Services2NObjMgr/EvalExp_Error: linking DlFunctionOfX"<<endl;
330 cout<<"...expval = "<<expval<<endl;
331 return;
332 }
333
334 // Evaluate function and close dynamic link
335 double result;
336 try {
337 result = f(0.);
338 CloseDLL();
339 } catch ( ... ) {
340 cout<<"Services2NObjMgr/EvalExp_Error: Arithmetic exception"<<endl;
341 CloseDLL();
342 return;
343 }
344
345 // Eventually store the result into variable or just print it
346 if(resultvarname.size()>0) {
347 if(MyObjMgr()->HasVar(resultvarname)) MyObjMgr()->DeleteVar(resultvarname);
348 char str[512];
349 if(result==0.) sprintf(str,"%f",result);
350 else {
351 double lr = log10(fabs(result));
352 if(lr>-17. && lr<17.) sprintf(str,"%17f",result);
353 else sprintf(str,"%.17e",result);
354 }
355 MyObjMgr()->SetVar(resultvarname,(string)str);
356 } else cout<<result<<" = "<<expval<<endl;
357
358}
359
360/* --Methode-- */
361void Services2NObjMgr::DisplayPoints2D(string& nom, string& expx, string& expy,
362 string& experrx, string& experry,
363 string& expcut, string dopt, string loop)
364{
365NObjMgrAdapter* obja=NULL;
366obja = MyObjMgr()->GetObjAdapter(nom);
367if (obja == NULL) {
368 cout << "Services2NObjMgr::DisplayPoints2D() Error , No such object " << nom << endl;
369 return;
370 }
371if (!mImgapp) return;
372
373// Creation NTuple
374const char* ntn[4] = {"expx","expy","expex","expey",};
375NTuple* nt = NULL;
376bool haserr = false;
377
378if ( (experrx.length() > 0 ) && (experry.length() > 0 ) ) { haserr = true; nt = new NTuple(4, ntn); }
379else { haserr = false; experrx = experry = "0."; nt = new NTuple(2, ntn); }
380
381ComputeExpressions(obja, expx, expy, experrx, experry, expcut, loop, nt, NULL, NULL);
382
383if (nt->NEntry() < 1) {
384 cout << "Services2NObjMgr::DisplayPoints2D() Warning Zero points satisfy cut !" << endl;
385 delete nt;
386 return;
387 }
388
389// nt->Show();
390// nt->Print(0,10);
391PINTuple* pin = new PINTuple(nt, true);
392pin->SelectXY(ntn[0], ntn[1]);
393if ( haserr ) pin->SelectErrBar(ntn[2], ntn[3]);
394
395dopt = "defline " + dopt;
396string titre = nom + ":" + expy + "%" + expx;
397mImgapp->DispScDrawer( (PIDrawer*)pin, titre, dopt);
398return;
399}
400
401/* --Methode-- */
402void Services2NObjMgr::DisplayPoints3D(string& nom, string& expx, string& expy,
403 string& expz, string& expcut, string dopt, string loop)
404{
405NObjMgrAdapter* obja=NULL;
406obja = MyObjMgr()->GetObjAdapter(nom);
407if (obja == NULL) {
408 cout << "Services2NObjMgr::DisplayPoints3D() Error , No such object " << nom << endl;
409 return;
410 }
411if (!mImgapp) return;
412
413const char* ntn[3] = {"expx","expy","expz"};
414NTuple* nt = new NTuple(3,ntn); // Creation NTuple
415
416string expwt = "1.";
417ComputeExpressions(obja, expx, expy, expz, expwt, expcut, loop, nt, NULL, NULL);
418
419if (nt->NEntry() < 1) {
420 cout << "Services2NObjMgr::DisplayPoints3D() Warning Zero points satisfy cut !" << endl;
421 delete nt;
422 return;
423 }
424//DBG nt->Show();
425//DBG nt->Print(0,10);
426PINTuple3D* pin = new PINTuple3D(nt, true);
427pin->SelectXYZ(ntn[0], ntn[1], ntn[2]);
428dopt = "defline " + dopt;
429
430// Pour plot a partir de DispScDrawer
431// string nomdisp = "_NT3D_";
432// mImgapp->DispScDrawer( (PIDrawer*)pin, nomdisp, opt);
433// Pour plot a partir de Disp3DDrawer
434string titre = nom + ":" + expy + "%" + expx;
435mImgapp->Disp3DDrawer(pin, titre, dopt);
436
437return;
438}
439
440/* --Methode-- */
441void Services2NObjMgr::DisplayPoints2DW(string& nom, string& expx, string& expy,
442 string& expwt, string& expcut, bool fgcolidx, string dopt, string loop)
443{
444NObjMgrAdapter* obja=NULL;
445obja = MyObjMgr()->GetObjAdapter(nom);
446if (obja == NULL) {
447 cout << "Services2NObjMgr::DisplayPoints2DW() Error , No such object " << nom << endl;
448 return;
449 }
450if (!mImgapp) return;
451
452const char* ntn[3] = {"expx","expy","expw"};
453NTuple* nt = new NTuple(3,ntn); // Creation NTuple
454
455string exp = "1.";
456ComputeExpressions(obja, expx, expy, expwt, exp, expcut, loop, nt, NULL, NULL);
457
458if (nt->NEntry() < 1) {
459 cout << "Services2NObjMgr::DisplayPoints2DW() Warning Zero points satisfy cut !" << endl;
460 delete nt;
461 return;
462 }
463
464PINTuple* pin = new PINTuple(nt, true);
465pin->SelectXY(ntn[0], ntn[1]);
466if (fgcolidx) pin->SelectColorByIndex(ntn[2]); // Use the weight expression as a color index
467else pin->SelectWt(ntn[2]);
468
469string titre = nom + ":" + expwt + "_" + expy + "%" + expx ;
470mImgapp->DispScDrawer( (PIDrawer*)pin, titre, dopt);
471return;
472}
473
474/* --Methode-- */
475void Services2NObjMgr::DisplayPoints3DW(string& nom, string& expx, string& expy, string& expz,
476 string& expwt, string& expcut, string dopt, string loop)
477{
478NObjMgrAdapter* obja=NULL;
479obja = MyObjMgr()->GetObjAdapter(nom);
480if (obja == NULL) {
481 cout << "Services2NObjMgr::DisplayPoints3D() Error , No such object " << nom << endl;
482 return;
483 }
484if (!mImgapp) return;
485
486const char* ntn[4] = {"expx","expy","expz","expw"};
487NTuple* nt = new NTuple(4,ntn); // Creation NTuple
488
489ComputeExpressions(obja, expx, expy, expz, expwt, expcut, loop, nt, NULL, NULL);
490
491if (nt->NEntry() < 1) {
492 cout << "Services2NObjMgr::DisplayPoints3DW() Warning Zero points satisfy cut !" << endl;
493 delete nt;
494 return;
495 }
496//DBG nt->Show();
497//DBG nt->Print(0,10);
498PINTuple3D* pin = new PINTuple3D(nt, true);
499pin->SelectXYZ(ntn[0], ntn[1], ntn[2]);
500pin->SelectWt(ntn[3]);
501pin->UseColorScale(true);
502dopt = "defline " + dopt;
503
504// Pour plot a partir de DispScDrawer
505// string nomdisp = "_NT3D_";
506// mImgapp->DispScDrawer( (PIDrawer*)pin, nomdisp, opt);
507// Pour plot a partir de Disp3DDrawer
508string titre = nom + ":" + expy + "%" + expx + " W=" + expwt;
509mImgapp->Disp3DDrawer(pin, titre, dopt);
510
511return;
512}
513
514/* --Methode-- */
515void Services2NObjMgr::ProjectH1(string& nom, string& expx, string& expwt,
516 string& expcut, string& nomh1, string dopt, string loop)
517{
518NObjMgrAdapter* obja=NULL;
519obja = MyObjMgr()->GetObjAdapter(nom);
520if (obja == NULL) {
521 cout << "Services2NObjMgr::ProjectH1() Error , No such object " << nom << endl;
522 return;
523}
524if (!mImgapp) return;
525
526Histo* h1 = NULL;
527NTuple* nt = NULL;
528AnyDataObj* oh = NULL;
529bool h1_already_exist = false;
530if (nomh1.length() > 0) oh=MyObjMgr()->GetObj(nomh1);
531else nomh1 = "/tmp/projh1d";
532if ( (oh != NULL) && (typeid(*oh) == typeid(Histo)) ) {
533 h1 = (Histo*)oh;
534 h1_already_exist = true;
535 // Pas de remise a zero ! h1->Zero();
536} else {
537 const char* ntn[2]= {"hxval", "hwt"};
538 nt = new NTuple(2,ntn); // Creation NTuple
539}
540string expz = "0.";
541ComputeExpressions(obja, expx, expwt, expz, expwt, expcut, loop, nt, h1, NULL);
542
543if ((!h1) && (!nt)) return;
544if (!h1) {
545 if (nt->NEntry() < 1) {
546 cout << "Services2NObjMgr::ProjectH1() Warning Zero points satisfy cut !" << endl;
547 delete nt;
548 return;
549 }
550 double xmin, xmax;
551 nt->GetMinMax(0, xmin, xmax);
552 h1 = new Histo(xmin, xmax, 100);
553 int k;
554 float* xn;
555 for(k=0; k<nt->NEntry(); k++) {
556 xn = nt->GetVec(k);
557 h1->Add(xn[0], xn[1]);
558 }
559 delete nt;
560 MyObjMgr()->AddObj(h1, nomh1);
561}
562
563if(!h1_already_exist || dopt.size()>0) MyObjMgr()->DisplayObj(nomh1, dopt);
564return;
565}
566
567/* --Methode-- */
568void Services2NObjMgr::ProjectH2(string& nom, string& expx, string& expy, string& expwt,
569 string& expcut, string& nomh2, string dopt, string loop)
570{
571NObjMgrAdapter* obja=NULL;
572obja = MyObjMgr()->GetObjAdapter(nom);
573if (obja == NULL) {
574 cout << "Services2NObjMgr::ProjectH2() Error , No such object " << nom << endl;
575 return;
576}
577if (!mImgapp) return;
578
579Histo2D* h2 = NULL;
580NTuple* nt = NULL;
581AnyDataObj* oh = NULL;
582bool h2_already_exist = false;
583if (nomh2.length() > 0) oh=MyObjMgr()->GetObj(nomh2);
584else nomh2 = "/tmp/projh2d";
585if ( (oh != NULL) && (typeid(*oh) == typeid(Histo2D)) ) {
586 h2 = (Histo2D*)oh;
587 h2_already_exist = true;
588 // Pas de remise a zero ! h2->Zero();
589} else {
590 const char* ntn[3]= {"hxval", "hyval", "hwt"};
591 nt = new NTuple(3,ntn); // Creation NTuple
592}
593string expz = "0.";
594ComputeExpressions(obja, expx, expy, expwt, expwt, expcut, loop, nt, NULL, h2);
595
596if ((!h2) && (!nt)) return;
597if (!h2) {
598 if (nt->NEntry() < 1) {
599 cout << "Services2NObjMgr::ProjectH2() Warning Zero points satisfy cut !" << endl;
600 delete nt;
601 return;
602 }
603 double xmin, xmax, ymin, ymax;
604 nt->GetMinMax(0, xmin, xmax);
605 nt->GetMinMax(1, ymin, ymax);
606 h2 = new Histo2D(xmin, xmax, 50, ymin, ymax, 50);
607 int k;
608 float* xn;
609 for(k=0; k<nt->NEntry(); k++) {
610 xn = nt->GetVec(k);
611 h2->Add(xn[0], xn[1], xn[2]);
612 }
613 delete nt;
614 MyObjMgr()->AddObj(h2, nomh2);
615}
616
617if(!h2_already_exist || dopt.size()>0) MyObjMgr()->DisplayObj(nomh2, dopt);
618return;
619
620}
621
622/* --Methode-- cmv 13/10/98 */
623void Services2NObjMgr::ProjectHProf(string& nom, string& expx, string& expy, string& expwt,
624 string& expcut, string& nomprof, string dopt, string loop)
625// Pour remplir un ``GeneralFitData'' a partir de divers objets:
626//| nom = nom de l'objet a projeter dans un HProf.
627//| expx = expression X de definition du bin.
628//| expy = expression Y a additionner dans le bin.
629//| expwt = expression W du poids a additionner.
630//| expcut = expression du test de selection.
631//| nomprof = nom du HProf engendre (optionnel). Si l'objet n'existe pas
632//| les limites Xmin,Xmax sont calculees automatiquement.
633//| sinon ce sont celles de l'objet preexistant.
634//| opt = options generales pour le display.
635{
636NObjMgrAdapter* obja=NULL;
637obja = MyObjMgr()->GetObjAdapter(nom);
638if (obja == NULL) {
639 cout << "Services2NObjMgr::ProjectHProf() Error , No such object " << nom << endl;
640 return;
641}
642if (!mImgapp) return;
643
644HProf* hprof = NULL;
645NTuple* nt = NULL;
646AnyDataObj* oh = NULL;
647bool hp_already_exist = false;
648if (nomprof.length() > 0) oh=MyObjMgr()->GetObj(nomprof);
649else nomprof = "/tmp/projprof";
650if( (oh!=NULL) && (typeid(*oh) == typeid(HProf)) ) {
651 hprof = (HProf*)oh;
652 hp_already_exist = true;
653} else {
654 const char* ntn[3]= {"hxval", "hyval", "hwt"};
655 nt = new NTuple(3,ntn); // Creation NTuple
656}
657string expz = "0.";
658ComputeExpressions(obja, expx, expy, expwt, expwt, expcut, loop, nt, NULL, NULL, hprof);
659
660if((!hprof) && (!nt)) return;
661if(!hprof) {
662 if (nt->NEntry() < 1) {
663 cout << "Services2NObjMgr::ProjectHProf() Warning Zero points satisfy cut !" << endl;
664 delete nt;
665 return;
666 }
667 r_8 xmin, xmax;
668 nt->GetMinMax(0, xmin, xmax);
669 hprof = new HProf(xmin, xmax, 100);
670 int k;
671 float* xn;
672 for(k=0; k<nt->NEntry(); k++) {
673 xn = nt->GetVec(k);
674 hprof->Add(xn[0], xn[1], xn[2]);
675 }
676 delete nt;
677 MyObjMgr()->AddObj(hprof, nomprof);
678}
679hprof->UpdateHisto();
680
681if(!hp_already_exist || dopt.size()>0) MyObjMgr()->DisplayObj(nomprof, dopt);
682return;
683}
684
685
686/* --Methode-- */
687void Services2NObjMgr::FillVect(string& nom, string& expx, string& expv,
688 string& expcut, string& nomvec, string dopt, string loop)
689{
690NObjMgrAdapter* obja=NULL;
691obja = MyObjMgr()->GetObjAdapter(nom);
692if (obja == NULL) {
693 cout << "Services2NObjMgr::FillVect() Error , No such object: " << nom << endl;
694 return;
695 }
696if (!mImgapp) return;
697
698Vector* v1 = NULL;
699AnyDataObj* ov = NULL;
700ov=MyObjMgr()->GetObj(nomvec);
701if (ov != NULL) v1 = dynamic_cast<Vector *>(ov);
702if (v1 == NULL) {
703 cout << "Services2NObjMgr::FillVect() Error , No such object or not a vector: " << nomvec << endl;
704 return;
705 }
706
707const char* ntn[2]= {"vi", "vv"};
708NTuple* nt = new NTuple(2,ntn); // Creation NTuple
709
710string expz = "0.";
711ComputeExpressions(obja, expx, expv, expz, expz, expcut, loop, nt);
712
713if (!nt) return;
714if (nt->NEntry() < 1) {
715 cout << "Services2NObjMgr::FillVect() Warning Zero points satisfy cut !" << endl;
716 delete nt;
717 return;
718 }
719
720 int i,k;
721 double* xn;
722 for(k=0; k<nt->NEntry(); k++) {
723 xn = nt->GetLineD(k);
724 i = int(xn[0]+0.5);
725 if ( (i < 0) || i >= v1->NElts() ) continue;
726 (*v1)(i) = xn[1];
727 }
728 delete nt;
729
730
731MyObjMgr()->DisplayObj(nomvec, dopt);
732return;
733}
734
735/* --Methode-- */
736void Services2NObjMgr::FillMatx(string& nom, string& expx, string& expy, string& expv,
737 string& expcut, string& nommtx, string dopt, string loop)
738{
739NObjMgrAdapter* obja=NULL;
740obja = MyObjMgr()->GetObjAdapter(nom);
741if (obja == NULL) {
742 cout << "Services2NObjMgr::FillMatx() Error , No such objet " << nom << endl;
743 return;
744 }
745if (!mImgapp) return;
746
747Matrix* mtx = NULL;
748AnyDataObj* om = NULL;
749om=MyObjMgr()->GetObj(nommtx);
750if (om != NULL) mtx = dynamic_cast<Matrix *>(om);
751if (mtx == NULL) {
752 cout << "Services2NObjMgr::FillMatx() Error , No such object or not a matrix " << nommtx << endl;
753 return;
754 }
755
756const char* ntn[3]= {"mi", "mj", "mv"};
757NTuple* nt = new NTuple(3,ntn); // Creation NTuple
758
759string expz = "0.";
760ComputeExpressions(obja, expx, expy, expv, expz, expcut, loop, nt);
761
762if (!nt) return;
763if (nt->NEntry() < 1) {
764 cout << "Services2NObjMgr::FillMatx() Warning Zero points satisfy cut !" << endl;
765 delete nt;
766 return;
767 }
768
769 int ic, jl, k;
770 double* xn;
771 for(k=0; k<nt->NEntry(); k++) {
772 xn = nt->GetLineD(k);
773 ic = int(xn[0]+0.5);
774 jl = int(xn[1]+0.5);
775 if ( (ic < 0) || ic >= mtx->NCol() ) continue;
776 if ( (jl < 0) || jl >= mtx->NRows() ) continue;
777 (*mtx)(jl, ic) = xn[2];
778 }
779 delete nt;
780
781
782MyObjMgr()->DisplayObj(nommtx, dopt);
783return;
784
785}
786
787/* --Methode-- */
788void Services2NObjMgr::ExpressionToVector(string& nom, string& expx, string& expcut,
789 string& nomvec, string dopt, string loop)
790{
791NObjMgrAdapter* obja=NULL;
792obja = MyObjMgr()->GetObjAdapter(nom);
793if (obja == NULL) {
794 cout << "Services2NObjMgr::ExpressionToVector() Error , No such object " << nom << endl;
795 return;
796 }
797if (!mImgapp) return;
798
799NTuple* nt = NULL;
800// if (nomvec.length() < 1) nomvec = "/tmp/expvec";
801
802const char* ntn[2]= {"vecval", "vecwt"};
803nt = new NTuple(1,ntn); // Creation NTuple
804
805string expwt = "1.";
806string expz = "0.";
807string dumexpcut = expcut; if(dumexpcut.size()<=0) dumexpcut = "1.";
808ComputeExpressions(obja,expx,expz,expz,expwt,dumexpcut,loop,nt,NULL,NULL);
809
810if (!nt) return;
811if (nt->NEntry() < 1) {
812 cout << "Services2NObjMgr::ExpressionToVector() Warning Zero points satisfy cut !" << endl;
813 delete nt;
814 return;
815 }
816
817Vector* vec = new Vector(nt->NEntry());
818int k;
819float* xn;
820for(k=0; k<nt->NEntry(); k++) {
821 xn = nt->GetVec(k);
822 (*vec)(k) = xn[0];
823 }
824delete nt;
825if (nomvec.size() > 0) {
826 MyObjMgr()->AddObj(vec, nomvec);
827 MyObjMgr()->DisplayObj(nomvec, dopt);
828}
829else {
830// On calcule et on affiche mean/sigma + min/max
831 double min, max, mean, sigma;
832 vec->MinMax(min, max);
833 MathArray<r_8> ma;
834 ma.MeanSigma(*vec, mean, sigma);
835 cout << "Size= " << vec->Size() << " Mean= " << mean << " Sigma= " << sigma
836 << " Min= " << min << " Max= " << max << endl;
837 delete vec;
838}
839return;
840}
841
842/* --Methode-- */
843void Services2NObjMgr::NtFromASCIIFile(string& nom,string& filename,double def_val)
844// Pour remplir un NTuple/DataTable "nom" existant a partir du fichier
845// ASCII table "filename". Si il y a plus de variables dans le
846// ntuple que dans le fichier "filename",
847// les sur-numeraires numeriques sont mises a "def_val" par defaut.
848{
849AnyDataObj* mobj = MyObjMgr()->GetObj(nom);
850if(mobj == NULL) {
851 cout<<"NtFromASCIIFile() Error, object "<<nom<<" not existing"<<endl;
852 return;
853}
854
855NTuple* nt = NULL;
856DataTable* dt = NULL;
857if(typeid(*mobj) == typeid(NTuple)) nt = (NTuple*) mobj;
858else if(typeid(*mobj) == typeid(DataTable)) dt = (DataTable*) mobj;
859
860if(nt==NULL && dt==NULL) {
861 cout<<"NtFromASCIIFile() Error, object "<<nom<<" not an NTuple nor a DataTable"<<endl;
862 return;
863}
864if (!mImgapp) return;
865
866if(nt) nt->FillFromASCIIFile(filename, def_val);
867else if(dt) {
868 ifstream is(filename.c_str());
869 dt->FillFromASCIIFile(is);
870}
871return;
872}
873
874/* --Methode-- */
875void Services2NObjMgr::FillNT(string& nom, string& expx, string& expy, string& expz,
876 string& expt, string& expcut, string& nomnt, string loop)
877{
878NObjMgrAdapter* obja=NULL;
879obja = MyObjMgr()->GetObjAdapter(nom);
880if (obja == NULL) {
881 cout << "Services2NObjMgr::FillNT() Error , No such object " << nom << endl;
882 return;
883 }
884if (!mImgapp) return;
885
886bool fgnnt = false;
887NTuple* nt = NULL;
888AnyDataObj* oh = NULL;
889if (nomnt.length() > 0) oh=MyObjMgr()->GetObj(nomnt);
890else nomnt = "/tmp/fillnt";
891if ( (oh != NULL) && (typeid(*oh) == typeid(NTuple)) ) {
892 nt = (NTuple*)oh;
893 if (nt->NVar() > 10) {
894 cout << "Services2NObjMgr::FillNT() Warning , Max 10 var in NTuple -> new NTuple" << endl;
895 nt = NULL;
896 }
897 }
898if (nt == NULL) {
899 const char* ntn[4]= {"x", "y","z","t"};
900 nt = new NTuple(4,ntn); // Creation NTuple
901 fgnnt = true;
902 }
903
904ComputeExpressions(obja, expx, expy, expz, expt, expcut, loop, nt, NULL, NULL);
905
906if (fgnnt) MyObjMgr()->AddObj(nt, nomnt);
907return;
908
909}
910
911/* --Methode-- */
912void Services2NObjMgr::FillNTFrCFile(string & nom, string const & fname,
913 string const & funcname, string & nomnt, string loop)
914{
915if (!mImgapp) return;
916
917NObjMgrAdapter* obja=NULL;
918obja = MyObjMgr()->GetObjAdapter(nom);
919if (obja == NULL) {
920 cout << "Services2NObjMgr::FillNTFrCFile( " << nom << "...) No such object" <<endl;
921 return;
922 }
923bool adel = true;
924NTupleInterface* objnt = obja->GetNTupleInterface(adel);
925if (objnt == NULL) {
926 cout << "Services2NObjMgr::FillNTFrCFile( " << nom << "...) Not an NTupleInterface !" <<endl;
927 return;
928 }
929
930// Pour synchronisation d'execution simultanee
931ZSync zs(mutx_dynlink); zs.NOp();
932
933NTLoopExprFunc f = (NTLoopExprFunc)LinkFunctionFromFile(fname, funcname);
934if (!f) {
935 cerr << "Services2NObjMgr::FillNTFrCFile Error Creation NTLoopExprFunc" << endl;
936 if (adel) delete objnt; // Delete de l'objet NTupleInterface si necessaire
937 return;
938 }
939
940bool fgnnt = false;
941NTuple* nt = NULL;
942if (nomnt.length() > 0) {
943 AnyDataObj* oh = NULL;
944 oh=MyObjMgr()->GetObj(nomnt);
945 if ( (oh != NULL) && (typeid(*oh) == typeid(NTuple)) ) {
946 nt = (NTuple*)oh;
947 if (nt->NVar() > 10) {
948 cout << "Services2NObjMgr::FillNTFrCFile() Warning , Max 10 var in NTuple -> new NTuple" << endl;
949 nt = NULL;
950 }
951 }
952 if (nt == NULL) {
953 const char* ntn[4]= {"x", "y","z","t"};
954 nt = new NTuple(4,ntn); // Creation NTuple
955 fgnnt = true;
956 }
957 }
958
959double xnt[10];
960
961int i,k;
962for(i=0; i<10; i++) xnt[i] = 0.;
963
964int_8 k1,k2,dk;
965k1 = 0; k2 = objnt->NbLines(); dk = 1;
966DecodeLoopParameters(loop, k1, k2, dk);
967if (k1 < 0) k1 = 0;
968if (k2 < 0) k2 = objnt->NbLines();
969if (k2 > (int_8)objnt->NbLines()) k2 = objnt->NbLines();
970if (dk <= 0) dk = 1;
971
972try {
973 double* xn;
974 int_8 kstart = k1, kend = k2;
975 for(k=kstart; k<kend; k+=dk) {
976 xn = objnt->GetLineD(k);
977 if (f((int_8_exprf)k, xn, xnt, xnt+1, xnt+2, xnt+3, (int_8_exprf)kstart,(int_8_exprf) kend) != 0) {
978 if (nt) nt->Fill(xnt);
979 }
980 }
981 }
982catch ( PException exc ) {
983 fflush(stdout);
984 cout << endl;
985 cerr << endl;
986 cerr << "Services2NObjMgr::FillNTFrCFile() Exception :" << exc.Msg() << endl;
987}
988
989if (adel) delete objnt; // Delete de l'objet NTupleInterface si necessaire
990CloseDLL();
991
992if (fgnnt) MyObjMgr()->AddObj(nt, nomnt);
993return;
994}
995
996/* --Methode-- */
997void Services2NObjMgr::PrepareNTExpressionCFile(string & nom, string const & fname,
998 string const & funcname)
999{
1000NObjMgrAdapter* obja=NULL;
1001obja = MyObjMgr()->GetObjAdapter(nom);
1002if (obja == NULL) {
1003 cout << "Services2NObjMgr::PrepareNTExpressionCFile( " << nom << "...) No such object" <<endl;
1004 return;
1005 }
1006bool adel = true;
1007NTupleInterface* objnt = obja->GetNTupleInterface(adel);
1008if (objnt == NULL) {
1009 cout << "Services2NObjMgr::PrepareNTExpressionCFile( " << nom
1010 << "...) No NTupleInterface !" <<endl;
1011 return;
1012 }
1013string vardec = objnt->VarList_C("_xnti_");
1014
1015FILE *fip;
1016if ((fip = fopen(fname.c_str(), "w")) == NULL) {
1017 cout << "Services2NObjMgr::PrepareNTExpressionCFile()_Error: fopen " << fname << endl;
1018 if (adel) delete objnt; // Delete de l'objet NTupleInterface si necessaire
1019 return;
1020 }
1021
1022// constitution du fichier des decalarations des variables de l'interface NTuple
1023fputs("#include <stdlib.h> \n", fip);
1024fputs("#include <stdio.h> \n", fip);
1025fputs("#include <math.h> \n\n", fip);
1026
1027fputs("/* ------ Compare bits on double --------- */ \n", fip);
1028fputs("typedef long long int_8_exprf;\n", fip);
1029fputs("int_8_exprf BitCmp64(double v,int_8_exprf flg)\n", fip);
1030fputs("{return ((int_8_exprf)((v<0.) ? v-0.1 : v+0.1))&flg;}\n", fip);
1031fputs("/* ------ Some random number generators --------- */ \n", fip);
1032fputs("#if defined(__ppc__) && defined(__MACH__) \n",fip);
1033fputs("#include <limits.h> \n", fip);
1034fputs("#define drand48() ((double)(random())/LONG_MAX) \n",fip);
1035fputs("#endif \n",fip);
1036fputs("#define frand01() ( (float) drand48() ) \n", fip);
1037fputs("#define drand01() drand48() \n", fip);
1038fputs("#define rand01() drand48() \n", fip);
1039fputs("#define frandpm1() ( 2. * frand01() - 1.) \n", fip);
1040fputs("#define drandpm1() ( 2. * drand01() - 1.) \n", fip);
1041fputs("#define randpm1() ( 2. * drand01() - 1.) \n", fip);
1042fputs("double NorRand(void) \n", fip);
1043fputs(" { \n double x,A,B; \n LAB10: \n A = drand01(); \n", fip);
1044fputs(" if ( A == 0. ) goto LAB10; \n B = drand01(); \n", fip);
1045fputs(" x = sqrt(-2.*log(A))*cos(2.*M_PI*B); \n", fip);
1046fputs(" return(x); \n } \n", fip);
1047fputs("#define GauRand() NorRand() \n", fip);
1048fputs("#define gaurand() NorRand() \n\n", fip);
1049
1050fputs("/* Angle conversions */\n", fip);
1051fputs("double deg2rad(double x) { return (x*M_PI/180.); } \n", fip);
1052fputs("double rad2deg(double x) { return (x*180./M_PI); } \n", fip);
1053fputs("/* Theta(0..Pi, phi(0..2Pi) conversion to Molleweide X,Y */\n", fip);
1054fputs("double tetphi2mollX(double teta, double phi) { \n", fip);
1055fputs(" if (phi<=M_PI) return(180.-sin(acos(2*teta/M_PI-1.))*((M_PI-phi)*180./M_PI)); \n", fip);
1056fputs(" else return(180.+sin(acos(2.*teta/M_PI-1.))*((phi-M_PI)*180./M_PI)); } \n", fip);
1057fputs("double tetphi2mollY(double teta) { \n", fip);
1058fputs(" return (90.-(teta*180./M_PI)); } \n\n", fip);
1059fputs("/* Longitude(0..360), Latitude(-90..90) conversion to Molleweide X,Y */\n", fip);
1060fputs("double longlat2mollX(double lng, double lat) { \n", fip);
1061fputs(" if (lng<=180.) return(180.-sin(acos(lat/90.))*(180.-lng)); \n", fip);
1062fputs(" else return(180.+sin(acos(lat/90.))*(lng-180.)); } \n", fip);
1063fputs("double longlat2mollY(double lat) { return lat; } \n\n\n", fip);
1064
1065fputs("/* NTupleInterface Variable declaration - Generated by piapp */\n", fip);
1066fputs("/* -- Services2NObjMgr::PrepareNTExpressionCFile() -- */ \n", fip);
1067fputs("/* _nl line number or sequential index : _nstart <= _nl < _nend */ \n\n", fip);
1068fprintf(fip,"int %s(int_8_exprf _nl, double* _xnti_, double* _rx_, double* _ry_, double* _rz_, \n",
1069 funcname.c_str());
1070fprintf(fip," double* _rt_, int_8_exprf _nstart, int_8_exprf _nend) \n");
1071fprintf(fip, "{ \n %s \n", vardec.c_str());
1072fputs(" if (!1) { /* Cut Expression failed */ \n", fip);
1073fputs(" *_rx_ = *_ry_ = *_rz_ = *_rt_ = 0.; return(0);", fip);
1074fputs(" } \n /* Cut expression satisfied */ \n", fip);
1075fputs(" *_rx_ = 1.; \n *_ry_ = 1.; \n *_rz_ = 1.; \n *_rt_ = 1.; \n", fip);
1076fputs(" return(1); \n} \n", fip);
1077
1078fclose(fip);
1079
1080if (adel) delete objnt; // Delete de l'objet NTupleInterface si necessaire
1081return;
1082}
1083
1084/* --Methode-- cmv 13/10/98 */
1085void Services2NObjMgr::FillGFD(string& nom, string& expx, string& expy, string& expz,
1086 string& experr, string& expcut, string& nomgfd, string loop)
1087// Pour remplir un ``GeneralFitData'' a partir de divers objets:
1088//| nom = nom de l'objet a transcrire selon 1D: Z=f(X) ou 2D: Z=f(X,Y) .
1089//| Vector,Matrix,Histo,HProf,Histo2D,Image<T>,StarList,NTuple,GeneralFitData
1090//| expx = expression X du GeneralFitData (1er abscisse)
1091//| expy = expression Y du GeneralFitData (2sd abscisse si non "", Z=f(X,Y))
1092//| expz = expression Z du GeneralFitData (valeur de l'ordonnee)
1093//| experr = expression de l'erreur sur l'ordonnee Z
1094//| expcut = expression du test de selection
1095//| nomgfd = nom du GeneralFitData engendre (optionnel)
1096{
1097NObjMgrAdapter* obja=NULL;
1098obja = MyObjMgr()->GetObjAdapter(nom);
1099if (obja == NULL) {
1100 cout << "Services2NObjMgr::FillGFD() Error , No such object "<<nom<<endl;
1101 return;
1102 }
1103if(!mImgapp) return;
1104
1105// 2D ou 3D?
1106int nvar = 2;
1107if(expy.length()<=0) {nvar = 1; expy = "0.";}
1108
1109// Creation NTuple Buffer
1110const char* ntn[4]= {"x","y","f","e"};
1111NTuple*nt = new NTuple(4,ntn);
1112
1113// Remplissage NTuple buffer
1114ComputeExpressions(obja, expx, expy, expz, experr, expcut, loop, nt, NULL, NULL);
1115if(nt->NEntry() < 1)
1116 {cout<<"Services2NObjMgr::FillGFD() Warning Zero points satisfy cut !"<<endl;
1117 delete nt; return;}
1118
1119//Remplissage de la structure GeneraFitData
1120if (nt->NEntry() <= 0) {
1121 cout<<"Services2NObjMgr::FillGFD() Warning - NData= " << nt->NEntry() << endl;
1122 delete nt;
1123 return;
1124 }
1125
1126GeneralFitData* gfd = new GeneralFitData(nvar,nt->NEntry(),0);
1127int k;
1128float* xn;
1129for(k=0; k<nt->NEntry(); k++) {
1130 xn = nt->GetVec(k);
1131 gfd->AddData(xn,xn[2],xn[3]);
1132}
1133
1134// Menage et table d'objets
1135delete nt;
1136MyObjMgr()->AddObj(gfd, nomgfd);
1137return;
1138}
1139
1140/* --Methode-- cmv 12/07/00 */
1141void Services2NObjMgr::FillGFDfrVec(string nomgfd,string namx,string namy,string namz,string name)
1142// Pour remplir un ``GeneralFitData'' a partir de vecteurs
1143//| gdfrvec nomgd X Y ! !
1144//| gdfrvec nomgd X Y ! EY
1145//| gdfrvec nomgd X Y Z !
1146//| gdfrvec nomgd X Y Z EZ
1147//| - nomgfd = nom du generaldata a remplir
1148//| - namx = nom du vecteur contenant les valeurs X
1149//| - namy = nom du vecteur contenant les valeurs Y
1150//| - namz = nom du vecteur contenant les valeurs Z (ou "!")
1151//| - name = nom du vecteur contenant les valeurs des erreurs EY ou EZ
1152{
1153// Decodage des noms des vecteurs pour le remplissage du generaldata
1154if(nomgfd=="!" || nomgfd.length()<1)
1155 {cout<<"FillGFDfrVec_Error: bad GenaralData name "<<nomgfd<<endl; return;}
1156if(namx=="!" || namx.length()<1)
1157 {cout<<"FillGFDfrVec_Error: bad X vector name "<<namx<<endl; return;}
1158if(namy=="!" || namy.length()<1)
1159 {cout<<"FillGFDfrVec_Error: bad Y vector name "<<namy<<endl; return;}
1160if(namz.length()<1) namz = "!";
1161if(name.length()<1) name = "!";
1162int nvar = 0;
1163if(namz=="!") nvar = 1; else nvar = 2;
1164
1165// Identify data
1166NamedObjMgr omg;
1167AnyDataObj* mobj = NULL;
1168Vector* v;
1169int nel = 0;
1170r_8 *x=NULL, *y=NULL, *z=NULL,* ez=NULL;
1171
1172if( (mobj=omg.GetObj(namx)) == NULL) {
1173 cout<<"FillGFDfrVec_Error: unknown X object "<<namx<<endl; return;
1174} else {
1175 v = (Vector*) mobj; x = v->Data(); nel=v->NElts();
1176}
1177
1178if( (mobj=omg.GetObj(namy)) == NULL) {
1179 cout<<"FillGFDfrVec_Error: unknown Y object "<<namy<<endl; return;
1180} else {
1181 v = (Vector*) mobj; y = z = v->Data(); if(v->NElts()<nel) nel=v->NElts();
1182}
1183
1184if( nvar==2 && (mobj=omg.GetObj(namz)) == NULL) {
1185 cout<<"FillGFDfrVec_Error: unknown Z object "<<namz<<endl; return;
1186} else {
1187 v = (Vector*) mobj; z = v->Data(); if(v->NElts()<nel) nel=v->NElts();
1188}
1189
1190if(name!="!") {
1191 if( (mobj=omg.GetObj(name)) == NULL) {
1192 cout<<"FillGFDfrVec_Error: unknown EZ object "<<name<<endl; return;
1193 } else {
1194 v = (Vector*) mobj; ez = v->Data(); if(v->NElts()<nel) nel=v->NElts();
1195 }
1196}
1197
1198if(nel<=0)
1199 {cout<<"FillGFDfrVec_Error: bad number of elements "<<nel<<endl; return;}
1200
1201// Create GeneralData and fill it with vectors
1202GeneralFitData* gfd = new GeneralFitData(nvar,nel+5,0);
1203if(nvar==1) gfd->SetData1(nel,x,z,ez); // On remplit Y=f(X)
1204else gfd->SetData2(nel,x,y,z,ez); // On remplit Z=f(X,y)
1205
1206// Menage et table d'objets
1207if( omg.GetObj(nomgfd) != NULL ) omg.DelObj(nomgfd);
1208MyObjMgr()->AddObj(gfd,nomgfd);
1209return;
1210}
1211
1212/* --Methode-- */
1213void Services2NObjMgr::ComputeExpressions(NObjMgrAdapter* obja, string& expx,
1214 string& expy, string& expz, string& expt, string& expcut, string& loop,
1215 NTuple* nt, Histo* h1, Histo2D* h2, HProf* hp)
1216{
1217if (obja == NULL) return;
1218bool adel = true;
1219NTupleInterface* objnt = obja->GetNTupleInterface(adel);
1220if (objnt == NULL) return;
1221string vardec = objnt->VarList_C("_zz6qi_");
1222
1223// Pour synchronisation d'execution simultanee
1224ZSync zs(mutx_dynlink); zs.NOp();
1225
1226PlotExprFunc f = LinkExprFunc(vardec, expx, expy, expz, expt, expcut);
1227if (!f) {
1228 cerr << "Services2NObjMgr::::ComputeExpressions() Error Creation PlotExprFunc " << endl;
1229 if (adel) delete objnt; // Delete de l'objet NTupleInterface si necessaire
1230 return;
1231 }
1232
1233double xnt[10];
1234
1235int_8 k;
1236for(k=0; k<10; k++) xnt[k] = 0.;
1237int_8 k1,k2,dk;
1238k1 = 0; k2 = objnt->NbLines(); dk = 1;
1239DecodeLoopParameters(loop, k1, k2, dk);
1240if (k1 < 0) k1 = 0;
1241if (k2 < 0) k2 = objnt->NbLines();
1242if (k2 > (int_8)objnt->NbLines()) k2 = objnt->NbLines();
1243if (dk <= 0) dk = 1;
1244
1245try {
1246 double* xn;
1247 for(k=k1; k<k2; k += dk) {
1248 xn = objnt->GetLineD(k);
1249 if (f((int_8_exprf)k,xn, xnt, xnt+1, xnt+2, xnt+3) != 0) {
1250 if (nt) nt->Fill(xnt);
1251 if (h1) h1->Add(xnt[0], xnt[3]);
1252 if (h2) h2->Add(xnt[0], xnt[1], xnt[3]);
1253 if (hp) hp->Add(xnt[0], xnt[1], xnt[3]);
1254 }
1255 }
1256 }
1257catch ( PException exc ) {
1258 fflush(stdout);
1259 cout << endl;
1260 cerr << endl;
1261 cerr << "Services2NObjMgr::ComputeExpressions() Exception :" << exc.Msg() << endl;
1262}
1263
1264if (adel) delete objnt; // Delete de l'objet NTupleInterface si necessaire
1265// Fermeture du fichier .so
1266CloseDLL();
1267return;
1268}
1269
1270
1271/* --Methode-- */
1272PlotExprFunc Services2NObjMgr::LinkExprFunc(string& vardec, string& expx, string& expy,
1273 string& expz, string& expt, string& cut)
1274{
1275FILE *fip;
1276string fname = TmpDir + "expf_pia_dl.c";
1277string cmd;
1278int rc;
1279
1280cmd = "rm -f " + fname;
1281rc = system(cmd.c_str());
1282//DBG printf("LinkExprFunc_Do> %s (Rc=%d)\n", cmd.c_str(), rc);
1283
1284if ((fip = fopen(fname.c_str(), "w")) == NULL) {
1285 string sn = fname;
1286 cout << "Services2NObjMgr/LinkExprFunc_Error: fopen( " << sn << endl;
1287 return(NULL);
1288 }
1289
1290// constitution du fichier a compiler
1291fputs("#include <stdlib.h> \n", fip);
1292fputs("#include <math.h> \n", fip);
1293
1294fputs("/* ------ Compare bits on double --------- */ \n", fip);
1295fputs("typedef long long int_8_exprf;\n", fip);
1296fputs("int_8_exprf BitCmp64(double v,int_8_exprf flg)\n", fip);
1297fputs("{return ((int_8_exprf)((v<0.) ? v-0.1 : v+0.1))&flg;}\n", fip);
1298fputs("/* ------ Some random number generators --------- */ \n", fip);
1299fputs("#if defined(__ppc__) && defined(__MACH__) \n",fip);
1300fputs("#include <limits.h> \n", fip);
1301fputs("#define drand48() ((double)(random())/LONG_MAX) \n",fip);
1302fputs("#endif \n",fip);
1303fputs("#define frand01() ( (float) drand48() ) \n", fip);
1304fputs("#define drand01() drand48() \n", fip);
1305fputs("#define rand01() drand48() \n", fip);
1306fputs("#define frandpm1() ( 2. * frand01() - 1.) \n", fip);
1307fputs("#define drandpm1() ( 2. * drand01() - 1.) \n", fip);
1308fputs("#define randpm1() ( 2. * drand01() - 1.) \n", fip);
1309fputs("double NorRand(void) \n", fip);
1310fputs(" { \n double x,A,B; \n LAB10: \n A = drand01(); \n", fip);
1311fputs(" if ( A == 0. ) goto LAB10; \n B = drand01(); \n", fip);
1312fputs(" x = sqrt(-2.*log(A))*cos(2.*M_PI*B); \n", fip);
1313fputs(" return(x); \n } \n", fip);
1314fputs("#define GauRand() NorRand() \n", fip);
1315fputs("#define gaurand() NorRand() \n\n", fip);
1316
1317fputs("/* Angle conversions */\n", fip);
1318fputs("double deg2rad(double x) { return (x*M_PI/180.); } \n", fip);
1319fputs("double rad2deg(double x) { return (x*180./M_PI); } \n", fip);
1320fputs("/* Theta(0..Pi, phi(0..2Pi) conversion to Molleweide X,Y */\n", fip);
1321fputs("double tetphi2mollX(double teta, double phi) { \n", fip);
1322fputs(" if (phi<=M_PI) return(180.-sin(acos(2*teta/M_PI-1.))*((M_PI-phi)*180./M_PI)); \n", fip);
1323fputs(" else return(180.+sin(acos(2.*teta/M_PI-1.))*((phi-M_PI)*180./M_PI)); } \n", fip);
1324fputs("double tetphi2mollY(double teta) { \n", fip);
1325fputs(" return (90.-(teta*180./M_PI)); } \n\n", fip);
1326fputs("/* Longitude(0..360), Latitude(-90..90) conversion to Molleweide X,Y */\n", fip);
1327fputs("double longlat2mollX(double lng, double lat) { \n", fip);
1328fputs(" if (lng<=180.) return(180.-sin(acos(lat/90.))*(180.-lng)); \n", fip);
1329fputs(" else return(180.+sin(acos(lat/90.))*(lng-180.)); } \n", fip);
1330fputs("double longlat2mollY(double lat) { return lat; } \n\n\n", fip);
1331
1332fputs("int expf_pia_dl_func(int_8_exprf _nl, double* _zz6qi_, double* _rx_6q_, double* _ry_6q_, double* _rz_6q_, double* _rt_6q_) \n{\n", fip);
1333fprintf(fip,"%s \n", vardec.c_str());
1334fprintf(fip, "if (!(%s)) { *_rx_6q_ = *_ry_6q_ = *_rz_6q_ = *_rt_6q_ = 0.; return(0); } \n", cut.c_str());
1335fprintf(fip, "*_rx_6q_ = %s ; \n", expx.c_str());
1336fprintf(fip, "*_ry_6q_ = %s ; \n", expy.c_str());
1337fprintf(fip, "*_rz_6q_ = %s ; \n", expz.c_str());
1338fprintf(fip, "*_rt_6q_ = %s ; \n", expt.c_str());
1339fputs("return(1); \n} \n", fip);
1340fclose(fip);
1341string func = "expf_pia_dl_func";
1342return((PlotExprFunc)LinkFunctionFromFile(fname, func));
1343}
1344
1345
1346/* --Methode-- */
1347DlFunction Services2NObjMgr::LinkFunctionFromFile(string const & fname, string const & funcname)
1348{
1349// Le link dynamique
1350CloseDLL();
1351dynlink = PDynLinkMgr::BuildFromCFile(fname);
1352if (dynlink == NULL) {
1353 cerr << "Services2NObjMgr/LinkFunctionFromFile_Erreur: Erreur creation/Ouverture SO " << endl;
1354 return(NULL);
1355 }
1356
1357DlFunction retfunc = dynlink->GetFunction(funcname);
1358if (retfunc == NULL) {
1359 string sn = funcname;
1360 cerr << "Services2NObjMgr/LinkExprFunc_Erreur: Erreur linking " << sn << endl;
1361 CloseDLL();
1362 return(NULL);
1363 }
1364else return(retfunc);
1365}
1366
1367/* --Methode-- */
1368void Services2NObjMgr::CloseDLL()
1369{
1370if (dynlink) delete dynlink; dynlink = NULL;
1371}
1372
1373/* --Methode-- */
1374int Services2NObjMgr::ExecuteCommand(string line)
1375{
1376 if (mImgapp == NULL) return(99);
1377 return(mImgapp->CmdInterpreter()->Interpret(line));
1378}
1379
1380// Fonction static
1381/* --Methode-- */
1382void Services2NObjMgr::DecodeLoopParameters(string& loop, int_8& i1, int_8& i2, int_8& di)
1383{
1384// Decode des paramatres de boucle for(int i=i1; i<i2; i+=di) specifies
1385// sous forme i1[:i2[:di]]
1386// cout << "LoopParam() " << loop << " I1=" << i1 << " I2=" << i2 << " DI=" << di;
1387size_t l = loop.length();
1388if (l < 1) return;
1389size_t p = loop.find(':');
1390if (p >= l) { i1 = atol(loop.c_str()); return; }
1391i1 = atol(loop.substr(0, p).c_str());
1392string aa = loop.substr(p+1);
1393p = aa.find(':');
1394if (p < aa.length() ) {
1395 i2 = atol(aa.substr(0,p).c_str());
1396 di = atol(aa.substr(p+1).c_str());
1397 }
1398else i2 = atol(aa.c_str());
1399// cout << "-> I1= " << i1 << " I2= " << i2 << " DI= " << di << endl;
1400return;
1401}
1402
1403/* --Methode-- */
1404string Services2NObjMgr::FileName2Name(string const & fn)
1405{
1406
1407char fsep[2] = {FILESEP, '\0'};
1408char tsep[2] = {'.', '\0'};
1409size_t p = fn.find_last_of(fsep);
1410size_t l = fn.length();
1411if (p >= l) p = 0;
1412else p++;
1413size_t q = fn.find_first_of(tsep,p);
1414if (q < p) q = l;
1415return(fn.substr(p,q-p));
1416}
1417
1418
1419const char* Services2NObjMgr::PClassIdToClassName(int cid)
1420{
1421 return("AnyDataObj");
1422}
Note: See TracBrowser for help on using the repository browser.