source: Sophya/trunk/SophyaPI/PIext/pawexecut.cc@ 2792

Last change on this file since 2792 was 2792, checked in by cmv, 20 years ago

spiapp cmd "newdt" pour creer une DataTable
spiapp cmd "n/read" modifier pour remplir aussi des DataTable (option -dt)
NtFromASCIIFile modifier pour remplir aussi des DataTable

mais le decodage FillFromASCII de la classe DataTAble est encore a faire

cmv 01/06/2005

File size: 67.9 KB
Line 
1#include <stdio.h>
2#include <stdlib.h>
3#include <ctype.h>
4#include <iostream>
5#include <typeinfo>
6
7#include "sopnamsp.h"
8#include "strutil.h"
9#include "strutilxx.h"
10#include "histos.h"
11#include "histos2.h"
12#include "hisprof.h"
13#include "histerr.h"
14#include "ntuple.h"
15#include "datatable.h"
16
17#include "pawexecut.h"
18#include "nobjmgr.h"
19#include "servnobjm.h"
20#include "nomgadapter.h"
21#include "pistdimgapp.h"
22#include "pihisto.h"
23
24#ifdef SANS_EVOLPLANCK
25#include "cvector.h"
26#include "matrix.h"
27#else
28#include "tmatrix.h"
29#include "tvector.h"
30#endif
31
32/* Reza + cmv 13/10/99 */
33
34inline bool __vector_sort_up_r8(r_8 x1,r_8 x2) {return x1<x2;}
35inline bool __vector_sort_down_r8(r_8 x1,r_8 x2) {return x1>x2;}
36inline bool __vector_sort_up_r4(r_4 x1,r_4 x2) {return x1<x2;}
37inline bool __vector_sort_down_r4(r_4 x1,r_4 x2) {return x1>x2;}
38inline bool __vector_sort_up_i4(int_4 x1,int_4 x2) {return x1<x2;}
39inline bool __vector_sort_down_i4(int_4 x1,int_4 x2) {return x1>x2;}
40inline bool __vector_sort_up_u8(uint_8 x1,uint_8 x2) {return x1<x2;}
41inline bool __vector_sort_down_u8(uint_8 x1,uint_8 x2) {return x1>x2;}
42
43uint_4 PAWExecutor::autoc_counter_ = 0;
44
45/* methode */
46PAWExecutor::PAWExecutor(PIACmd *piac, PIStdImgApp* app)
47: mApp(app)
48{
49string kw, usage;
50string hgrp = "pawCmd";
51
52kw = "reset";
53usage = "Reset histograms vectors or matrix";
54usage += "\n reset nameobj";
55piac->RegisterCommand(kw,usage,this,hgrp);
56
57kw = "n/plot";
58usage = "Plot NTuple variables \"a la paw\" (alias n/pl)";
59usage += "\n n/plot nameobj.x_exp [cut] [w_exp] [loop] [gratt]";
60usage += "\n n/plot nameobj.y_exp%x_exp [cut] [loop] [gratt]";
61usage += "\n n/plot nameobj.z_exp%y_exp%x_exp [cut] [loop] [gratt]";
62usage += "\n for default use ! , loop=i1[:i2[:di]]";
63usage += "\n for 1 dimensional (1D) projection:";
64usage += "\n use graphic option \"keepbin\" to keep previous";
65usage += "\n plot binning for 1D distribution";
66usage += "\n Related commands: plot2dw plot3d";
67piac->RegisterCommand(kw,usage,this,hgrp);
68
69kw = "n/pl";
70usage = "alias to n/plot (see n/plot)";
71piac->RegisterCommand(kw,usage,this,hgrp);
72
73kw = "n/proj";
74usage = "Project NTuple in histogram (1D or 2D) a la paw";
75usage += "\n n/proj nameproj nameobj.x_exp [cut] [w_exp] [loop] [gratt]";
76usage += "\n n/proj nameproj nameobj.y_exp%x_exp [cut] [w_exp] [loop] [gratt]";
77usage += "\n for default use ! , loop=i1[:i2[:di]]";
78usage += "\n for 1 dimensional (1D) projection:";
79usage += "\n no display is performed if \"nameproj\" is an existing histogram";
80usage += "\n unless a graphic option \"gratt\" is given";
81usage += "\n Related commands: projh1d projh2d projprof exptovec";
82piac->RegisterCommand(kw,usage,this,hgrp);
83
84kw = "n/scan";
85usage = "Scan NTuple a la paw";
86usage += "\n n/scan nameobj[.exp1%exp2%exp3] cut loop";
87usage += "\n [-f:filename] [list_of_variables]";
88usage += "\n loop : iev1[:iev2[:diev]] or !";
89usage += "\n cut : cut expression or 1. or !";
90usage += "\n list_of_variables : default is all variables";
91usage += "\n : var1 var2 var3 ... varn";
92usage += "\n : var1 : var2 (from var1 to var2)";
93usage += "\n : : var2 (from first variable to var2)";
94usage += "\n : var1 : (from var1 to last variable)";
95usage += "\n ex: \"v1 : v3 v7 v4 : v6 v2 v9 :\"";
96usage += "\n exp1%exp2%exp3 :";
97usage += "\n if given add exp1,exp2,exp3 to the variable list";
98usage += "\n -f:filename : write into \"filename\", Default is to stdout";
99piac->RegisterCommand(kw,usage,this,hgrp);
100
101kw = "n/read";
102usage = "Read columns in an ASCII file and fill a NTuple or a DataTable";
103usage += "\n n/read nt fascii [options] var_1,c_1 var_2,c_2 ... var_n,c_n ";
104usage += "\n nt : NTuple or DataTable (see options) name";
105usage += "\n var_i,c_i : variable name, associated column in ASCII file [0,n[";
106usage += "\n where [options] are:";
107usage += "\n \"-dt\": create a DataTable of double instead of a NTuple";
108usage += "\n \"=s\": separator character is \'s\' (could be \"\t\")";
109usage += "\n \"-^abcd\": do not read lines beginning with string \"abcd\" ";
110usage += "\n \"+^abcd\": read only lines beginning with string \"abcd\" ";
111usage += "\n \"-abcd\": do not read lines which contain string \"abcd\" ";
112usage += "\n \"+abcd\": read only lines which contain string \"abcd\" ";
113usage += "\n these options may be repeated (ex: \"-^abcd\" \"-^xyz\") ";
114usage += "\n - in case of \"do not read\" options are added with logical AND ";
115usage += "\n - in case of \"read only\" options are added with logical OR ";
116piac->RegisterCommand(kw,usage,this,hgrp);
117
118kw = "n/merge";
119usage = "Merge ntuples";
120usage += "\n n/merge nt nt_1 nt_2 ... nt_n";
121usage += "\n Merge ntuples nt_i into ntuple nt";
122piac->RegisterCommand(kw,usage,this,hgrp);
123
124kw = "h/integ";
125usage = "Integrate a 1D histogram";
126usage += "\n h/integ nameh1d [norm]";
127usage += "\n norm<=0 means no normalisation (def norm=1)";
128usage += "\n Related commands: h/deriv v/integ v/deriv";
129piac->RegisterCommand(kw,usage,this,hgrp);
130
131kw = "v/integ";
132usage = "Integrate a TVector / vector";
133usage += "\n v/integ namevec [norm]";
134usage += "\n norm<=0 means no normalisation (def norm=0)";
135usage += "\n Related commands: h/integ h/deriv v/deriv";
136piac->RegisterCommand(kw,usage,this,hgrp);
137
138kw = "v/sort";
139usage = "Sort a vector into itself";
140usage += "\n v/sort namevec [+1/-1]";
141usage += "\n \"+1\" means increasing order, \"-1\" decreasing order";
142piac->RegisterCommand(kw,usage,this,hgrp);
143
144kw = "h/deriv";
145usage = "Derivate a 1D histogram";
146usage += "\n h/deriv nameh1d";
147usage += "\n Related commands: h/integ v/integ v/deriv";
148piac->RegisterCommand(kw,usage,this,hgrp);
149
150kw = "v/deriv";
151usage = "Derivate a TVector / vector";
152usage += "\n v/deriv namevec [deriv_option]";
153usage += "\n deriv_option -1 replace v[i] with v[i]-v[i-1]";
154usage += "\n 0 replace v[i] with (v[i+1]-v[i-1])/2 (default)";
155usage += "\n +1 replace v[i] with v[i+1]-v[i]";
156usage += "\n Related commands: h/integ h/deriv v/integ";
157piac->RegisterCommand(kw,usage,this,hgrp);
158
159kw = "h/rebin";
160usage = "Rebin a 1D histogram or profile";
161usage += "\n h/rebin nameh1d nbin";
162piac->RegisterCommand(kw,usage,this,hgrp);
163
164kw = "h/cadd";
165usage = "Add a constant to an histogram, a vector or a matrix";
166usage += "\n h/cadd namehisto val";
167usage += "\n Related commands: h/cmult h/oper";
168piac->RegisterCommand(kw,usage,this,hgrp);
169
170kw = "h/cmult";
171usage = "Multiply an histogram, a vector or a matrix by a constant";
172usage += "\n h/cmult namehisto val";
173usage += "\n Related commands: h/cadd h/oper";
174piac->RegisterCommand(kw,usage,this,hgrp);
175
176kw = "h/oper";
177usage = "Operation on histograms vectors or matrices";
178usage += "\n h/oper @ h1 h2 hres";
179usage += "\n hres = h1 @ h2 with @ = (+,-,*,/)";
180usage += "\n For vectors and matrices, operations are elements by elements";
181usage += "\n Related commands: h/cadd h/cmult";
182piac->RegisterCommand(kw,usage,this,hgrp);
183
184kw = "h/plot/2d";
185usage = "Specific plot for 2D histogrammes";
186usage += "\n h/plot/2d nameh2d show : infos on 2D histogramme";
187usage += "\n h/plot/2d nameh2d h [dopt] : plot 2D histogramme";
188usage += "\n h/plot/2d nameh2d px [dopt] : plot X projection";
189usage += "\n h/plot/2d nameh2d py [dopt] : plot Y projection";
190usage += "\n h/plot/2d nameh2d bx n [dopt] : plot X band number n";
191usage += "\n h/plot/2d nameh2d by n [dopt] : plot Y band number n";
192usage += "\n h/plot/2d nameh2d sx n [dopt] : plot X slice number n";
193usage += "\n h/plot/2d nameh2d sy n [dopt] : plot Y slice number n";
194usage += "\n n < 0 means Show Info";
195piac->RegisterCommand(kw,usage,this,hgrp);
196
197kw = "h/put_vec";
198usage = "Put content of vector (matrix) into content of histogram 1D or 2D";
199usage += "\n h/put_vec nameh1d namevector [cont,err2]";
200usage += "\n h/put_vec nameh2d namematrix [cont,err2]";
201usage += "\n cont : put into histogramme content";
202usage += "\n err2 : put into histogramme error^2";
203usage += "\n Related commands: h/get_vec";
204piac->RegisterCommand(kw,usage,this,hgrp);
205
206kw = "h/get_vec";
207usage = "Get content of histogram 1D profile or 2D into vector (matrix)";
208usage += "\n h/get_vec nameh1d namevector [cont,err2,absc] [proj]";
209usage += "\n h/get_vec nameh2d namematrix [cont,err2,absc]";
210usage += "\n cont : get histogramme content";
211usage += "\n err2 : get histogramme error^2";
212usage += "\n absc : get histogramme low bin abscissa (1D only)";
213usage += "\n proj :";
214usage += "\n show : show available projections for Histo2D";
215usage += "\n px : get X projection";
216usage += "\n py : get Y projection";
217usage += "\n bx n : get X band number n";
218usage += "\n by n : get Y band number n";
219usage += "\n sx n : get X slice number n";
220usage += "\n sy n : get Y slice number n";
221usage += "\n Related commands: h/put_vec";
222piac->RegisterCommand(kw,usage,this,hgrp);
223
224kw = "h/copy";
225usage = "Copy content of object1 into object2";
226usage += "\n objects are Vector,Matrix,Histo,Histo2D";
227usage += "\n h/copy namefrom nametocopy [i1[:i2]] [j1[:j2]] [ic1[:jc1]]";
228usage += "\n copy obj1Dfrom(i1->i2) into obj1Dto(ic1->)";
229usage += "\n copy obj2Dfrom(i1,j1->i2,j2) into obj2Dto(ic1,jc1->)";
230usage += "\n Warning: elements from i1 to i2 included are copied";
231usage += "\n If obj1Dto does not exist, is is created with size i2-i1+1";
232usage += "\n or obj2Dto with size i2-i1+1,j2-j1+1";
233usage += "\n The adressed content of obj?Dfrom is overwritten";
234usage += "\n The non-adressed content of obj?Dfrom is left unchanged";
235usage += "\n Related commands: copy";
236piac->RegisterCommand(kw,usage,this,hgrp);
237
238kw = "h/set/err";
239usage = "Set Histo,Histo2D errors for range of bins or range of values";
240usage += "\n h/set/err namehisto setvalue i1[:i2] [j1[:j2]]";
241usage += "\n set error to setvalue for bin range i1:i2 j1:j2";
242usage += "\n h/set/err namehisto setvalue v v1:v2";
243usage += "\n set error to setvalue for content values range v1:v2";
244usage += "\n h/set/err namehisto setvalue e e1:e2";
245usage += "\n set error to setvalue for error values range v1:v2";
246usage += "\n Related commands: h/set/cont";
247piac->RegisterCommand(kw,usage,this,hgrp);
248
249kw = "h/set/cont";
250usage = "Set Histo,Histo2D content for range of bins or range of values";
251usage += "\n h/set/cont namehisto setvalue i1[:i2] [j1[:j2]]";
252usage += "\n set content to setvalue for bin range i1:i2 j1:j2";
253usage += "\n h/set/cont namehisto setvalue v v1:v2";
254usage += "\n set content to setvalue for content values range v1:v2";
255usage += "\n h/set/cont namehisto setvalue e e1:e2";
256usage += "\n set content to setvalue for error values range v1:v2";
257usage += "\n Related commands: h/set/err";
258piac->RegisterCommand(kw,usage,this,hgrp);
259
260kw = "h/err";
261usage = "Set Histo,Histo2D error to function of bin content value";
262usage += "\n h/err namehisto expr_func";
263usage += "\n Related commands: h/set/err";
264piac->RegisterCommand(kw,usage,this,hgrp);
265
266}
267
268/* methode */
269PAWExecutor::~PAWExecutor()
270{
271}
272
273/* methode */
274int PAWExecutor::Execute(string& kw, vector<string>& tokens, string& toks)
275{
276if(kw == "reset") {
277 reset(tokens); return(0);
278} else if(kw == "n/plot" || kw == "n/pl") {
279 n_plot(tokens); return(0);
280} else if(kw == "n/proj") {
281 n_proj(tokens); return(0);
282} else if(kw == "n/scan") {
283 n_scan(tokens); return(0);
284} else if(kw == "n/read") {
285 n_read(tokens); return(0);
286} else if(kw == "n/merge") {
287 n_merge(tokens); return(0);
288} else if(kw == "h/integ") {
289 h_integ(tokens); return(0);
290} else if(kw == "v/integ") {
291 v_integ(tokens); return(0);
292} else if(kw == "v/sort") {
293 v_sort(tokens); return(0);
294} else if(kw == "h/deriv") {
295 h_deriv(tokens); return(0);
296} else if(kw == "v/deriv") {
297 v_deriv(tokens); return(0);
298} else if(kw == "h/rebin") {
299 h_rebin(tokens); return(0);
300} else if(kw == "h/cadd") {
301 h_cadd(tokens); return(0);
302} else if(kw == "h/cmult") {
303 h_cmult(tokens); return(0);
304} else if(kw == "h/oper") {
305 h_oper(tokens); return(0);
306} else if(kw == "h/plot/2d") {
307 h_plot_2d(tokens); return(0);
308} else if(kw == "h/put_vec") {
309 h_put_vec(tokens); return(0);
310} else if(kw == "h/get_vec") {
311 h_get_vec(tokens); return(0);
312} else if(kw == "h/copy") {
313 h_copy(tokens); return(0);
314} else if(kw == "h/set/err") {
315 string dum = "err";
316 h_set(dum,tokens); return(0);
317} else if(kw == "h/set/cont") {
318 string dum = "cont";
319 h_set(dum,tokens); return(0);
320} else if(kw == "h/err") {
321 h_err(tokens); return(0);
322} else return(1);
323}
324
325/* methode */
326bool PAWExecutor::IsThreadable(string const & keyw)
327{
328 if ( (keyw == "n/plot") || (keyw == "n/pl") ||
329 (keyw == "n/proj") ) return true;
330 else return false;
331}
332
333
334/* methode */
335void PAWExecutor::reset(vector<string>& tokens)
336// Reset d'histogrammes, vecteurs et matrices
337{
338if(tokens.size() < 1)
339 {cout<<"Usage: reset nameobj"<<endl; return;}
340NamedObjMgr omg;
341AnyDataObj* mobj = omg.GetObj(tokens[0]);
342if(mobj == NULL)
343 {cout<<"PAWExecutor::reset Error , Pas d'objet de nom "<<tokens[0]<<endl;
344 return;}
345string ctyp = typeid(*mobj).name();
346
347#ifdef SANS_EVOLPLANCK
348if(typeid(*mobj)==typeid(Vector)) {Vector* ob=(Vector*) mobj; ob->Zero();}
349else if(typeid(*mobj)==typeid(Matrix)) {Matrix* ob=(Matrix*) mobj; ob->Zero();}
350#else
351 if(typeid(*mobj)==typeid(Vector)) {Vector* ob=(Vector*) mobj; (*ob) = 0.; }
352// ob->DataBlock().Reset(0.);}
353 else if(typeid(*mobj)==typeid(Matrix)) {Matrix* ob=(Matrix*) mobj; (*ob) = 0.; }
354//ob->DataBlock().Reset(0.);}
355#endif
356else if(typeid(*mobj)==typeid(Histo)) {Histo* ob=(Histo*) mobj; ob->Zero();}
357else if(typeid(*mobj)==typeid(HProf)) {HProf* ob=(HProf*) mobj; ob->Zero();}
358else if(typeid(*mobj)==typeid(HistoErr)) {HistoErr* ob=(HistoErr*) mobj; ob->Zero();}
359else if(typeid(*mobj)==typeid(Histo2D)) {Histo2D* ob=(Histo2D*)mobj; ob->Zero();}
360else {
361 cout<<"PAWExecutor::reset Error , No reset possible on "<<ctyp<<endl;
362 return;
363}
364
365return;
366}
367
368/* methode */
369void PAWExecutor::n_plot(vector<string>& tokens)
370// Equivalent n/plot de paw
371// Plot 1D
372// n/plot nameobj.x_exp [cut] [w_exp] [gratt]
373// Plot 2D (plot2dw)
374// n/plot nameobj.y_exp%x_exp [cut] [w_exp] [gratt]
375// Plot 3D (plot3d)
376// n/plot nameobj.z_exp%y_exp%x_exp [cut] [gratt]
377{
378if(tokens.size() < 1) {
379 cout
380 <<"Usage: n/plot nameobj.x_exp [cut] [w_exp] [loop] [gratt] [nomh1]"<<endl
381 <<" n/plot nameobj.y_exp%x_exp [cut] [loop] [gratt]"<<endl
382 <<" n/plot nameobj.z_exp%y_exp%x_exp [cut] [loop] [gratt]"<<endl
383 <<" for default use ! , loop=i1[:i2[:di]]"<<endl;
384 return;
385}
386string nameobj,expx,expy,expz;
387int_4 nvar = decodepawstring(tokens[0],nameobj,expx,expy,expz);
388string expcut = "1"; string expwt = "1."; string loop = "";
389string dopt = ""; string nameproj="";
390if(tokens.size()>=2) expcut = tokens[1]; if(expcut=="!") expcut="1";
391
392NamedObjMgr omg;
393Services2NObjMgr* srvo = omg.GetServiceObj();
394
395if(nvar<=0) {
396 cout<<"PAWExecutor::n_plot Error: bad coding "<<tokens[0]<<endl;
397} else if(nvar==1) { // c'est un plot 1D
398 if(tokens.size()>=3) expwt = tokens[2]; if(expwt=="!") expwt="1.";
399 if(tokens.size()>=4) loop = tokens[3]; if(loop=="!") loop="";
400 if(tokens.size()>=5) dopt = tokens[4]; if(dopt=="!") dopt="";
401 if(tokens.size()>=6) nameproj = tokens[5];
402 if(nameproj.length()<=0 || nameproj=="!") {
403 nameproj = "/autoc/paw_n_plot1D_";
404 AnyDataObj* mobj = omg.GetObj(nameproj);
405 if(mobj!=NULL) {
406 char buff[16]; autoc_counter_++; sprintf(buff,"%d",autoc_counter_);
407 nameproj += buff;
408 }
409 // --- Si option "keepbin" on projete avec le meme binning que l'histo 1D precedent
410 if(dopt.find("keepbin")<dopt.size()) {
411 PIStdImgApp* piimapp = omg.GetImgApp();
412 if(piimapp) {
413 PIBaseWdg* pibwdg = piimapp->CurrentBaseWdg();
414 if(pibwdg) {
415 int nbdrw = pibwdg->NbDrawers();
416 if(nbdrw>0) {
417 for(int i=nbdrw-1;i>=0;i--) {
418 PIDrawer* pidwr = pibwdg->GetDrawer(i);
419 PIHisto* pih = NULL;
420 if( (pih = dynamic_cast<PIHisto *>(pidwr))==NULL ) continue;
421 Histo* h = pih->Histogram();
422 if(h==NULL) continue;
423 if(h->NBins()<1) continue;
424 Histo* hsame = new Histo(h->XMin(),h->XMax(),h->NBins());
425 omg.AddObj(hsame,nameproj);
426 break;
427 }
428 }
429 }
430 }
431 }
432 // ---
433 }
434 srvo->ProjectH1(nameobj,expx,expwt,expcut,nameproj,dopt,loop);
435} else if(nvar==2) { // c'est un plot 2D
436 if(tokens.size()>=3) loop = tokens[2]; if(loop=="!") loop="";
437 if(tokens.size()>=4) dopt = tokens[3];
438 string err = "";
439 srvo->DisplayPoints2D(nameobj,expx,expy,err,err,expcut,dopt,loop);
440} else { // c'est un plot 3D
441 if(tokens.size()>=3) loop = tokens[2]; if(loop=="!") loop="";
442 if(tokens.size()>=4) dopt = tokens[3];
443 srvo->DisplayPoints3D(nameobj,expx,expy,expz,expcut,dopt,loop);
444}
445
446return;
447}
448
449
450/* methode */
451void PAWExecutor::n_proj(vector<string>& tokens)
452// Equivalent n/proj de paw
453// Project NTuple in histogram a la paw
454// Dans un Histo 1D
455// n/proj nameproj nameobj.x_exp [cut] [w_exp] [gratt]
456// Dans un Histo 2D ou un HProf (dans ce cas nameproj doit etre cree).
457// n/proj nameproj nameobj.y_exp%x_exp [cut] [w_exp] [gratt]
458{
459if(tokens.size()<2)
460 {cout<<"Usage: n/proj nameproj nameobj.[y_exp%]x_exp [cut] [w_exp] [loop] [gratt]"<<endl
461 <<" for default use ! , loop=i1[:i2[:di]]"<<endl; return;}
462string nameproj = tokens[0];
463string nameobj,expx,expy,expz;
464int_4 nvar = decodepawstring(tokens[1],nameobj,expx,expy,expz);
465string expcut = "1"; string expwt = "1."; string loop = ""; string dopt = "";
466if(tokens.size()>=3) expcut = tokens[2]; if(expcut=="!") expcut="1";
467if(tokens.size()>=4) expwt = tokens[3]; if(expwt=="!") expwt="1.";
468if(tokens.size()>=5) loop = tokens[4]; if(loop=="!") loop="";
469if(tokens.size()>=6) dopt = tokens[5];
470
471NamedObjMgr omg;
472Services2NObjMgr* srvo = omg.GetServiceObj();
473
474if(nvar==1) {
475 // c'est une projection dans un histo 1D
476 srvo->ProjectH1(nameobj,expx,expwt,expcut,nameproj,dopt,loop);
477} else if(nvar==2) {
478 // c'est une projection dans un histo2D
479 // OU un HProf si nameproj est un HProf un deja defini
480 AnyDataObj* mobj = omg.GetObj(nameproj);
481 if(mobj==NULL)
482 srvo->ProjectH2(nameobj,expx,expy,expwt,expcut,nameproj,dopt,loop);
483 else if(dynamic_cast<HProf*>(mobj))
484 srvo->ProjectHProf(nameobj,expx,expy,expwt,expcut,nameproj,dopt,loop);
485 else
486 srvo->ProjectH2(nameobj,expx,expy,expwt,expcut,nameproj,dopt,loop);
487} else {
488 cout<<"PAWExecutor::n_proj Error: bad coding "<<tokens[1]<<" nvar="<<nvar<<endl;
489}
490
491return;
492}
493
494/* methode */
495void PAWExecutor::n_scan(vector<string>& tokens)
496{
497if(tokens.size()<3)
498 {cerr<<"Usage: n/scan nameobj[.exp1%exp2%exp3] cut loop\n"
499 <<" [-f:filename] [list_of_variables]"<<endl; return;}
500
501// decodage des premiers arguments
502string nameobj,expr[4];
503int_4 nexpr = decodepawstring(tokens[0],nameobj,expr[0],expr[1],expr[2]);
504if(nexpr<4) {for(int_4 i=nexpr;i<4;i++) expr[i]="1.";}
505string expcut = tokens[1]; if(expcut=="!") expcut="1";
506string loop = tokens[2]; if(loop=="!") loop="";
507FILE* fout = NULL;
508uint_4 ldebvar = 3;
509if(tokens.size()>3) {
510 if(tokens[3].find("-f:") == 0) {
511 string filename = tokens[3].substr(3);
512 if(filename.size()>0) {
513 fout = fopen(filename.c_str(),"w");
514 if(fout==NULL)
515 {cerr<<"PAWExecutor::n_scan Error, file "<<filename
516 <<" not opened"<<endl; return;}
517 }
518 ldebvar++;
519 }
520}
521
522// ntuple adaptateur
523NamedObjMgr omg;
524Services2NObjMgr& srvo = *omg.GetServiceObj();
525NObjMgrAdapter* obja = omg.GetObjAdapter(nameobj); // Ne pas deleter
526if(obja == NULL)
527 {cerr<<"PAWExecutor::n_scan Error, ObjAdapter==NULL for "
528 <<nameobj<<endl; return;}
529bool adel = true;
530NTupleInterface* objnt = obja->GetNTupleInterface(adel);
531if(objnt == NULL)
532 {cerr<<"PAWExecutor::n_scan Error, NTupleInterface==NULL for "
533 <<nameobj<<endl; return;}
534
535// function pour le choix
536string vardec = objnt->VarList_C("_zz6qi_");
537PlotExprFunc f = srvo.LinkExprFunc(vardec,expr[0],expr[1],expr[2],expr[3],expcut);
538if(!f)
539 {cerr<<"PAWExecutor::n_scan Error, Creation PlotExprFunc"<<endl;
540 if(adel) delete objnt; if(fout) fclose(fout); return;}
541
542// variables a imprimer
543// "rien" --> de 0 a ncol-1 (toutes les variables)
544// : v1 --> de 0 a v1
545// v1 : --> de v1 a ncol-1
546// v1 : v2 --> de v1 a v2 (si v2 apres v1)
547// v1 et v2 (si v2 avant v1)
548// v1 v2 v3 v4 --> v1 v2 v3 v4 (ordre indifferent)
549// et toute combinaison : "v1 : v3 v7 v4 : v6 v2 v9 :"
550// --> v1 v2 v3 v7 v4 v5 v6 v2 v9 v10...v(ncol-1)
551int_4 ncol = objnt->NbColumns();
552if(ncol<=0)
553 {cerr<<"PAWExecutor::n_scan Error, no columns for NTuple"<<endl;
554 return;}
555vector<int_4> varnum;
556if(ldebvar>=tokens.size()) { // Toutes les variables
557 for(int_4 i=0;i<ncol;i++) varnum.push_back(i);
558} else { // Choix de certaines variables
559 int_4 k,klast,kk; bool frlast=false;
560 if(tokens[ldebvar]==":") {varnum.push_back(0); frlast=true;}
561 else {k = objnt->ColumnIndex(tokens[ldebvar]); varnum.push_back(k);}
562 ldebvar++;
563 if(ldebvar<tokens.size()) for(uint_4 i=ldebvar;i<tokens.size();i++) {
564 if(tokens[i]!=":") { // pas un separateur
565 k = klast = objnt->ColumnIndex(tokens[i]);
566 if(frlast) klast = varnum[varnum.size()-1] + 1;
567 if(klast>k) klast=k;
568 for(kk=klast;kk<=k;kk++) varnum.push_back(kk);
569 frlast=false;
570 } else if(i==tokens.size()-1) { // separateur a la fin
571 k = ncol-1;
572 klast = varnum[varnum.size()-1] + 1;
573 if(klast>k) klast=k;
574 for(kk=klast;kk<=k;kk++) varnum.push_back(kk);
575 } else frlast=true; // separateur pas a la fin
576 }
577}
578
579vector<string> varname;
580if(varnum.size()>0) for(int_4 i=0;i<(int)varnum.size();i++) {
581 if(varnum[i]<0 || varnum[i]>=ncol)
582 {cerr<<"PAWExecutor::n_scan Error, bad variable name at pos "
583 <<i<<endl; if(adel) delete objnt; if(fout) fclose(fout); return;}
584 string dum = objnt->ColumnName(varnum[i]);
585 varname.push_back(dum);
586}
587
588// evenements a utiliser
589int_8 k1=0, k2=objnt->NbLines(), dk=1;
590srvo.DecodeLoopParameters(loop,k1,k2,dk);
591if (k1<0) k1=0;
592if (k2<0) k2=objnt->NbLines();
593if (k2>(int_8)objnt->NbLines()) k2=objnt->NbLines();
594if (dk<=0) dk=1;
595
596// boucle sur les evenements et print
597try {
598 int_4 i;
599 if(fout) fprintf(fout,"#ev "); else printf("#ev ");
600 for(i=0;i<(int)varname.size();i++)
601 if(fout) fprintf(fout,"%s ",varname[i].c_str());
602 else printf( "%s ",varname[i].c_str());
603 if(nexpr>0) for(i=0;i<nexpr;i++)
604 if(fout) fprintf(fout,"%s ",expr[i].c_str());
605 else printf( "%s ",expr[i].c_str());
606 if(fout) fprintf(fout,"\n"); else printf("\n");
607
608 double xnt[5]={0,0,0,0,0};
609 double* xn;
610 for(int_8 k=k1; k<k2; k += dk) {
611 xn = objnt->GetLineD(k);
612 if(f((int_8_exprf)k,xn,xnt,xnt+1,xnt+2,xnt+3) != 0) {
613 if(fout) fprintf(fout,"%ld ",(long)k); else printf("%ld ",(long)k);
614 for(i=0;i<(int)varnum.size();i++) {
615 if(fout) fprintf(fout,"%g ",*(xn+varnum[i]));
616 else printf( "%g ",*(xn+varnum[i]));
617 }
618 if(nexpr>0) for(i=0;i<nexpr;i++) {
619 if(fout) fprintf(fout,"%g ",*(xnt+i));
620 else printf( "%g ",*(xnt+i));
621 }
622 if(fout) fprintf(fout,"\n"); else printf("\n");
623 }
624 }
625} // fin du try
626#ifdef SANS_EVOLPLANCK
627CATCH(merr) {
628 fflush(stdout); cout<<endl; cerr<<endl;
629 string es = PeidaExc(merr);
630 cerr<<"Services2NObjMgr::ComputeExpressions() Exception :"<<merr<<es;
631} ENDTRY;
632#else
633catch ( PException exc ) {
634 fflush(stdout); cout<<endl; cerr<<endl;
635 cerr<<"Services2NObjMgr::ComputeExpressions() Exception :"<<exc.Msg()<<endl;
636}
637#endif
638
639if(adel) delete objnt;
640if(fout) fclose(fout);
641srvo.CloseDLL(); // Fermeture du fichier .so
642return;
643}
644
645
646/* methode */
647#define __LENLINE_N_READ__ 8192
648void PAWExecutor::n_read(vector<string>& tokens)
649{
650 int lp=1;
651
652 if(tokens.size()<3) {
653 cerr<<"Usage: n/read nt fascii [options] var_1,c_1 ... var_n,c_n"<<endl;
654 return;
655 }
656
657 // decodage des arguments
658 string nament = tokens[0];
659 string nameascii = tokens[1];
660 bool create_datatable = false;
661 vector<string> donotreadbeg;
662 vector<string> donotreadin;
663 vector<string> onlyreadbeg;
664 vector<string> onlyreadin;
665 vector<string> varname;
666 vector<int> colnum;
667 char separator = ' ';
668 int numcolmaxi=-1;
669
670 for(int i=2;i<tokens.size();i++) {
671 int lc = tokens[i].size();
672 if(lc<2) continue;
673 if(tokens[i].find("-dt")==0) { // create DataTable instead of NTuple
674 create_datatable = true;
675 continue;
676 }
677 const char *c = tokens[i].c_str();
678 if(c[0]=='=') { // Separator
679 separator = c[1];
680 if(lc==3) if(c[1]=='\\' && c[2]=='t') separator = '\t';
681 continue;
682 }
683 if(c[0]=='+') { // Selection des lignes a lire
684 if(c[1]!='^') onlyreadin.push_back(&c[1]);
685 else if(lc>2) onlyreadbeg.push_back(&c[2]);
686 continue;
687 }
688 if(c[0]=='-') { // Selection des lignes commentaire
689 if(c[1]!='^') donotreadin.push_back(&c[1]);
690 else if(lc>2) donotreadbeg.push_back(&c[2]);
691 continue;
692 }
693 // decodage des noms de variables et des colonnes associees
694 int p = tokens[i].find(',');
695 if(p<1 || p>=lc-1) continue;
696 string vn = tokens[i].substr(0,p);
697 string cn = tokens[i].substr(p+1,lc-p-1);
698 int ic = atoi(cn.c_str());
699 if(ic<0) continue;
700 if( !isalpha(vn[0]) ) continue;
701 if(ic>numcolmaxi) numcolmaxi = ic;
702 varname.push_back(vn);
703 colnum.push_back(ic);
704 }
705
706 int nvar = varname.size();
707 if(nvar<=0) {
708 cerr<<"n_read: no variables to be read"<<endl;
709 return;
710 }
711
712 // Print what has to be done
713 if(lp) {
714 if(onlyreadin.size()>0) {
715 cout<<"n_read Only read line containing ["<<onlyreadin.size()<<"]:";
716 for(int i=0;i<onlyreadin.size();i++) cout<<" \'"<<onlyreadin[i]<<"\'";
717 cout<<endl;
718 }
719 if(onlyreadbeg.size()>0) {
720 cout<<"n_read Only read line begining with ["<<onlyreadbeg.size()<<"]:";
721 for(int i=0;i<onlyreadbeg.size();i++) cout<<" \'"<<onlyreadbeg[i]<<"\'";
722 cout<<endl;
723 }
724 if(donotreadin.size()>0) {
725 cout<<"n_read Do not read line containing ["<<donotreadin.size()<<"]:";
726 for(int i=0;i<donotreadin.size();i++) cout<<" \'"<<donotreadin[i]<<"\'";
727 cout<<endl;
728 }
729 if(donotreadbeg.size()>0) {
730 cout<<"n_read Do not read line begining with ["<<donotreadbeg.size()<<"]:";
731 for(int i=0;i<donotreadbeg.size();i++) cout<<" \'"<<donotreadbeg[i]<<"\'";
732 cout<<endl;
733 }
734 if(nvar>0) {
735 cout<<"n_read Number of variables to be read: "<<nvar<<endl;
736 for(int i=0;i<nvar;i++) cout<<" \'"<<varname[i].c_str()<<","<<colnum[i]<<"\'";
737 cout<<endl;
738 }
739 if(separator!=' ') cout<<"n_read Separator is: \'"<<separator<<"\'"<<endl;
740 }
741
742 // Open ASCII file
743 FILE * fascii = fopen(nameascii.c_str(),"r");
744 if(fascii==NULL) {
745 cerr<<"n_read: cannot open file "<<nameascii<<endl;
746 return;
747 }
748
749 // Creation du NTuple ou de la DataTable
750 NTuple *nt = NULL;
751 DataTable *dt = NULL;
752 if(!create_datatable) { //create NTuple
753 char** ntvn = new char*[nvar];
754 for(int i=0;i<nvar;i++) ntvn[i] = const_cast<char *>(varname[i].c_str());
755 nt = new NTuple(nvar,ntvn);
756 delete [] ntvn;
757 } else { //create DataTable
758 dt = new DataTable();
759 for(int i=0;i<nvar;i++) dt->AddDoubleColumn(varname[i]);
760 }
761 r_8 *xnt = new r_8[nvar];
762
763 // Read file
764 char *line = new char[__LENLINE_N_READ__];
765 int nline=0, nlinecom=0;
766 while(fgets(line,__LENLINE_N_READ__,fascii) != NULL ) {
767 nline++;
768 int lc = strlen(line); if(lc<1) continue;
769 // Pour enlever le \n final
770 if(line[lc-1]=='\n' || line[lc-1]=='\r')
771 {line[lc-1]='\0'; lc = strlen(line); if(lc<1) continue;}
772
773 string const sline(line);
774 //cout<<"\'"<<sline<<"\' lc="<<lc<<endl;
775
776 // Faut t'il lire cette ligne ?
777 bool read_line_1 = true;
778 if(onlyreadin.size()>0 || onlyreadbeg.size()>0) read_line_1 = false;
779 if(onlyreadin.size()>0) {
780 for(int i=0;i<onlyreadin.size();i++) {
781 uint_4 p = sline.find(onlyreadin[i].c_str());
782 if(p>=lc) continue;
783 read_line_1 = true;
784 break;
785 }
786 }
787 if(onlyreadbeg.size()>0) {
788 for(int i=0;i<onlyreadbeg.size();i++) {
789 uint_4 p = sline.find(onlyreadbeg[i].c_str());
790 if(p!=0) continue;
791 read_line_1 = true;
792 break;
793 }
794 }
795
796 // Faut t'il ne pas lire cette ligne ?
797 bool read_line_2 = true;
798 if(donotreadin.size()>0) {
799 for(int i=0;i<donotreadin.size();i++) {
800 uint_4 p = sline.find(donotreadin[i].c_str());
801 if(p>=lc) continue;
802 read_line_2 = false;
803 break;
804 }
805 }
806 if(donotreadbeg.size()>0) {
807 for(int i=0;i<donotreadbeg.size();i++) {
808 uint_4 p = sline.find(donotreadbeg[i].c_str());
809 if(p!=0) continue;
810 read_line_2 = false;
811 break;
812 }
813 }
814 if(!read_line_2) nlinecom++;
815
816 if(!read_line_1 || !read_line_2) continue;
817
818 // Decodage de la ligne
819 vector<string> vs;
820 //FillVStringFrString(sline,vs,separator);
821 SplitStringToVString(sline,vs,separator);
822 int lvs = vs.size();
823 //for(int i=0;i<lvs;i++) cout<<"|"<<vs[i]<<"| "; cout<<endl;
824 if(lvs<numcolmaxi) continue; // Pas assez de champs decodes, mauvaise ligne
825
826 // Remplissage du NTuple ou de la DataTable
827 for(int i=0;i<nvar;i++) {
828 xnt[i] = 0.;
829 int ic = colnum[i];
830 if(ic>=lvs) continue;
831 xnt[i] = atof(vs[ic].c_str());
832 }
833 if(nt) nt->Fill(xnt);
834 else if(dt) dt->Fill(xnt);
835 //cout<<"...xnt"; for(int i=0;i<nvar;i++) cout<<" "<<xnt[i]; cout<<endl;
836
837 }
838 int_4 nentries = (nt)? nt->NEntry(): dt->NEntry();
839 cout<<"n_read: "<<nline<<" lines in file, "
840 <<nlinecom<<" commentary, "
841 <<nentries<<" object entries"<<endl;
842
843 // On sauve le NTuple/DataTable si besoin, on ferme et detruit ce qu'il faut
844 NamedObjMgr omg;
845 if(nentries>0) {
846 if(nt) omg.AddObj(nt,nament);
847 if(dt) omg.AddObj(dt,nament);
848 } else {
849 if(nt) delete nt;
850 if(dt) delete dt;
851 }
852 delete [] xnt;
853 delete [] line;
854 fclose(fascii);
855
856 return;
857}
858#undef __LENLINE_N_READ__
859
860/* methode */
861void PAWExecutor::n_merge(vector<string>& tokens)
862{
863 if(tokens.size()<2) {
864 cerr<<"Usage: n/read nt nt_1 nt_2 ... nt_n"<<endl;
865 return;
866 }
867
868 NamedObjMgr omg;
869
870 // decodage des arguments
871 string nament = tokens[0];
872
873 // boucle sur les ntuples
874 NTuple * nt = NULL;
875 int nvar=0, nfill=0;
876 r_8 *xnt=NULL;
877
878 for(int i=1;i<tokens.size();i++) {
879
880 AnyDataObj* mobj = omg.GetObj(tokens[i]);
881 if(mobj==NULL) {
882 cout<<"n_merge Error: unknow object"<<tokens[i]<<endl;
883 continue;;
884 }
885 NTuple* nt1 = dynamic_cast<NTuple*>(mobj);
886 if(nt1==NULL) {
887 cout<<"n_merge Error: "<<tokens[i]<<" not a NTuple"<<endl;
888 continue;
889 }
890 if(nt1->NEntry()==0) {
891 cout<<"n_merge Error: "<<tokens[i]<<" is empty"<<endl;
892 continue;
893 }
894 if(nt1->NVar()==0) {
895 cout<<"n_merge Error: "<<tokens[i]<<" has no variable"<<endl;
896 continue;
897 }
898
899 // create receiving ntuple if first pass
900 if(nt==NULL) {
901 nvar = nt1->NVar();
902 vector<string> sntvn;
903 for(int i=0;i<nvar;i++) sntvn.push_back(nt1->ColumnName(i));
904 char **ntvn = new char*[nvar];
905 for(int i=0;i<nvar;i++) ntvn[i] = const_cast<char *>(sntvn[i].c_str());
906 nt = new NTuple(nvar,ntvn);
907 delete [] ntvn;
908 }
909
910 // filling with current ntuple
911 int nvar1 = nt1->NVar();
912 int n = (nvar1>nvar)? nvar1: nvar;
913 r_8 *xnt1 = new r_8[n];
914 for(int i=0;i<n;i++) xnt1[i]=0.;
915 for(uint_4 iev=0;iev<nt1->NEntry();iev++) {
916 nt1->GetVecD(iev,xnt1);
917 nt->Fill(xnt1);
918 }
919 nfill++;
920 delete [] xnt1;
921
922 }
923
924 if(xnt!=NULL) delete [] xnt;
925 if(nt!=NULL) {
926 cout<<"n_merge: ntuple filled with "<<nfill
927 <<" ntuples, "<<nt->NEntry()<<" entries"<<endl;
928 if(nt->NEntry()>0) omg.AddObj(nt,nament); else delete nt;
929 }
930 return;
931}
932
933/* methode */
934void PAWExecutor::h_integ(vector<string>& tokens)
935// Pour remplacer le contenu d'un histo 1D par son integrale
936{
937if(tokens.size()<1)
938 {cout<<"Usage: h/integ nameh1d [norm]"<<endl; return;}
939NamedObjMgr omg;
940AnyDataObj* mobj = omg.GetObj(tokens[0]);
941if(mobj==NULL)
942 {cout<<"PAWExecutor::h_integ Error: unknow object"<<tokens[0]<<endl;
943 return;}
944r_8 norm = 1.;
945if(tokens.size()>=2) norm = atof(tokens[1].c_str());
946// attention: dynamic_cast<Histo*>(HProf)=Vrai!
947Histo* h1 = dynamic_cast<Histo*>(mobj);
948HProf* hp = dynamic_cast<HProf*>(mobj);
949if(hp || !h1)
950 {cout<<"PAWExecutor::h_integ Error: "<<tokens[0]<<" not an Histo"<<endl;
951 return;}
952h1->HInteg(norm);
953}
954
955/* methode */
956void PAWExecutor::v_integ(vector<string>& tokens)
957// Pour remplacer le contenu d'un TVector<r_8> par son integrale
958// normalisee a norm:
959// Si norm <= 0 : pas de normalisation, integration seule
960{
961if(tokens.size()<1)
962 {cout<<"Usage: v/integ namevec [norm]"<<endl; return;}
963NamedObjMgr omg;
964AnyDataObj* mobj = omg.GetObj(tokens[0]);
965if(mobj==NULL)
966 {cout<<"PAWExecutor::v_integ Error: unknow object"<<tokens[0]<<endl;
967 return;}
968r_8 norm=-1.;
969if(tokens.size()>=2) norm = atof(tokens[1].c_str());
970#ifdef SANS_EVOLPLANCK
971Vector* v = dynamic_cast<Vector*>(mobj);
972#else
973TVector<r_8>* v = dynamic_cast<TVector<r_8>*>(mobj);
974#endif
975
976if(!v)
977 {cout<<"PAWExecutor::v_integ Error: "<<tokens[0]<<" not a TVector/Vector"<<endl;
978 return;}
979
980#ifdef SANS_EVOLPLANCK
981uint_4 n = v->NElts();
982#else
983uint_4 n = v->Size();
984#endif
985
986if(n==0)
987 {cout<<"PAWExecutor::v_integ Error: "<<tokens[0]<<" is an empty vector"<<endl;
988 return;}
989if(n>1) for(uint_4 i=1;i<n;i++) (*v)(i)+=(*v)(i-1);
990if(norm<=0. || (*v)(n-1)==0.) return;
991norm /= (*v)(n-1);
992for(uint_4 i=0;i<n;i++) (*v)(i) *= norm;
993}
994
995/* methode */
996void PAWExecutor::v_sort(vector<string>& tokens)
997// Pour sort in-place d'un vecteur ascendant ou descendant
998{
999 if(tokens.size()<1)
1000 {cout<<"Usage: v/sort namevec [+1/-1]"<<endl; return;}
1001 NamedObjMgr omg;
1002 AnyDataObj* mobj = omg.GetObj(tokens[0]);
1003 if(mobj==NULL)
1004 {cout<<"PAWExecutor::v_sort Error: unknow object"<<tokens[0]<<endl;
1005 return;}
1006
1007 bool up = true;
1008 if(tokens.size()>=2) if(atoi(tokens[1].c_str())<0) up=false;
1009
1010 {
1011 TVector<r_8>* v = dynamic_cast<TVector<r_8>*>(mobj);
1012 if(v) {
1013 if(v->Size()==0) return;
1014 if(up) stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_up_r8);
1015 else stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_down_r8);
1016 cout<<"median is: "<<0.5*((*v)(v->Size()/2)+(*v)((v->Size()-1)/2))<<endl;
1017 return;
1018 }
1019 }
1020
1021 {
1022 TVector<r_4>* v = dynamic_cast<TVector<r_4>*>(mobj);
1023 if(v) {
1024 if(v->Size()==0) return;
1025 if(up) stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_up_r4);
1026 else stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_down_r4);
1027 cout<<"median is: "<<0.5*((*v)(v->Size()/2)+(*v)((v->Size()-1)/2))<<endl;
1028 return;
1029 }
1030 }
1031
1032 {
1033 TVector<int_4>* v = dynamic_cast<TVector<int_4>*>(mobj);
1034 if(v) {
1035 if(v->Size()==0) return;
1036 if(up) stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_up_i4);
1037 else stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_down_i4);
1038 cout<<"median is: "<<0.5*((*v)(v->Size()/2)+(*v)((v->Size()-1)/2))<<endl;
1039 return;
1040 }
1041 }
1042
1043 {
1044 TVector<uint_8>* v = dynamic_cast<TVector<uint_8>*>(mobj);
1045 if(v) {
1046 if(v->Size()==0) return;
1047 if(up) stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_up_u8);
1048 else stable_sort(v->Data(), v->Data()+v->Size(), __vector_sort_down_u8);
1049 cout<<"median is: "<<0.5*((*v)(v->Size()/2)+(*v)((v->Size()-1)/2))<<endl;
1050 return;
1051 }
1052 }
1053
1054 cout<<"PAWExecutor::v_sort Error: "<<tokens[0]
1055 <<" not a TVector or not a supported TVector<TYPE>"<<endl;
1056 return;
1057
1058}
1059
1060/* methode */
1061void PAWExecutor::h_deriv(vector<string>& tokens)
1062// Pour remplacer le contenu d'un histo 1D par sa derivee
1063{
1064if(tokens.size()<1)
1065 {cout<<"Usage: h/deriv nameh1d"<<endl; return;}
1066NamedObjMgr omg;
1067AnyDataObj* mobj = omg.GetObj(tokens[0]);
1068if(mobj==NULL)
1069 {cout<<"PAWExecutor::h_deriv Error: unknow object"<<tokens[0]<<endl;
1070 return;}
1071// attention: dynamic_cast<Histo*>(HProf)=Vrai!
1072Histo* h1 = dynamic_cast<Histo*>(mobj);
1073HProf* hp = dynamic_cast<HProf*>(mobj);
1074if(hp || !h1)
1075 {cout<<"PAWExecutor::h_deriv Error: "<<tokens[0]<<" not an Histo"<<endl;
1076 return;}
1077h1->HDeriv();
1078}
1079
1080/* methode */
1081void PAWExecutor::v_deriv(vector<string>& tokens)
1082// Pour remplacer le contenu d'un TVector<r_8> par sa derivee
1083// deriv_option = -1 replace v[i] with v[i]-v[i-1]
1084// = 0 replace v[i] with (v[i+1]-v[i-1])/2 (default)
1085// = +1 replace v[i] with v[i+1]-v[i]
1086{
1087if(tokens.size()<1)
1088 {cout<<"Usage: v/deriv namevec [deriv_option(-1,0,+1)]"<<endl; return;}
1089NamedObjMgr omg;
1090AnyDataObj* mobj = omg.GetObj(tokens[0]);
1091if(mobj==NULL)
1092 {cout<<"PAWExecutor::v_deriv Error: unknow object"<<tokens[0]<<endl;
1093 return;}
1094int_4 deriv_option = 0;
1095if(tokens.size()>=2) deriv_option = atoi(tokens[1].c_str());
1096
1097#ifdef SANS_EVOLPLANCK
1098Vector* v = dynamic_cast<Vector*>(mobj);
1099#else
1100TVector<r_8>* v = dynamic_cast<TVector<r_8>*>(mobj);
1101#endif
1102
1103if(!v)
1104 {cout<<"PAWExecutor::v_deriv Error: "<<tokens[0]<<" not a TVector/Vector"<<endl;
1105 return;}
1106
1107
1108#ifdef SANS_EVOLPLANCK
1109uint_4 n = v->NElts();
1110#else
1111uint_4 n = v->Size();
1112#endif
1113
1114if(n==0)
1115 {cout<<"PAWExecutor::v_deriv Error: "<<tokens[0]<<" is an empty vector"<<endl;
1116 return;}
1117if(n<=1) return;
1118
1119#ifdef SANS_EVOLPLANCK
1120Vector vsave(*v);
1121#else
1122TVector<r_8> vsave(*v,false);
1123#endif
1124
1125if(deriv_option<0) {
1126 for(uint_4 i=1;i<n;i++) (*v)(i) = vsave(i)-vsave(i-1);
1127 (*v)(0) = (*v)(1);
1128} else if(deriv_option>0) {
1129 for(uint_4 i=0;i<n-1;i++) (*v)(i) = vsave(i+1)-vsave(i);
1130 (*v)(n-1) = (*v)(n-2);
1131} else {
1132 for(uint_4 i=1;i<n-1;i++) (*v)(i) = (vsave(i+1)-vsave(i-1))/2.;
1133 (*v)(0) = vsave(1)-vsave(0);
1134 (*v)(n-1) = vsave(n-1)-vsave(n-2);
1135}
1136}
1137
1138
1139/* methode */
1140void PAWExecutor::h_rebin(vector<string>& tokens)
1141// Pour re-binner un histogramme 1D
1142{
1143if(tokens.size()<2)
1144 {cout<<"Usage: h/rebin nameh1d nbin"<<endl; return;}
1145NamedObjMgr omg;
1146AnyDataObj* mobj = omg.GetObj(tokens[0]);
1147if(mobj==NULL)
1148 {cout<<"PAWExecutor::h_rebin Error: unknow object"<<tokens[0]<<endl;
1149 return;}
1150int_4 nbin = atoi(tokens[1].c_str());
1151Histo* h1 = dynamic_cast<Histo*>(mobj);
1152// Ca marche aussi pour les HProf, HRebin a ete passe virtuel
1153//HProf* hp = dynamic_cast<HProf*>(mobj); if(hp || !h1)
1154if(!h1)
1155 {cout<<"PAWExecutor::h_rebin Error: "<<tokens[0]<<" not an Histo/HProf"<<endl;
1156 return;}
1157h1->HRebin(nbin);
1158}
1159
1160/* methode */
1161void PAWExecutor::h_cadd(vector<string>& tokens)
1162// Additionne une constante a un histogramme
1163{
1164if(tokens.size()<2)
1165 {cout<<"Usage: h/cadd nameh1d val"<<endl; return;}
1166NamedObjMgr omg;
1167AnyDataObj* mobj = omg.GetObj(tokens[0]);
1168if(mobj==NULL)
1169 {cout<<"PAWExecutor::h_cadd Error: unknow object"<<tokens[0]<<endl;
1170 return;}
1171r_8 val = atof(tokens[1].c_str());
1172Histo* h1 = dynamic_cast<Histo*>(mobj);
1173HProf* hp = dynamic_cast<HProf*>(mobj);
1174Histo2D* h2 = dynamic_cast<Histo2D*>(mobj);
1175Vector* v = dynamic_cast<Vector*>(mobj);
1176Matrix* m = dynamic_cast<Matrix*>(mobj);
1177if(h1 && !hp) *h1 += val;
1178else if(h2) *h2 += val;
1179else if(v) *v += val;
1180else if(m) *m += val;
1181else cout<<"PAWExecutor::h_cadd Error: not implemented for "<<typeid(*mobj).name()<<endl;
1182}
1183
1184/* methode */
1185void PAWExecutor::h_cmult(vector<string>& tokens)
1186// Multiplie un histogramme par une constante
1187{
1188if(tokens.size()<2)
1189 {cout<<"Usage: h/cmult nameh1d val"<<endl; return;}
1190NamedObjMgr omg;
1191AnyDataObj* mobj = omg.GetObj(tokens[0]);
1192if(mobj==NULL)
1193 {cout<<"PAWExecutor::h_cmult Error: unknow object"<<tokens[0]<<endl;
1194 return;}
1195r_8 val = atof(tokens[1].c_str());
1196Histo* h1 = dynamic_cast<Histo*>(mobj);
1197HProf* hp = dynamic_cast<HProf*>(mobj);
1198Histo2D* h2 = dynamic_cast<Histo2D*>(mobj);
1199Vector* v = dynamic_cast<Vector*>(mobj);
1200Matrix* m = dynamic_cast<Matrix*>(mobj);
1201if(h1 && !hp) *h1 *= val;
1202else if(h2) *h2 *= val;
1203else if(v) *v += val;
1204else if(m) *m += val;
1205else cout<<"PAWExecutor::h_mult Error: not implemented for "<<typeid(*mobj).name()<<endl;
1206}
1207
1208/* methode */
1209void PAWExecutor::h_oper(vector<string>& tokens)
1210// Equivalent h/oper/add sub,mul,div de paw
1211// Operation entre 2 histogrammes ou 2 vecteurs ou 2 matrices
1212// h/oper @ h1 h2 hres
1213// hres = h1 @ h2 with @ = (+,-,*,/)
1214// Pour les vecteurs et les matrices, il sagit d'operations elements a elements
1215{
1216if(tokens.size()<4)
1217 {cout<<"Usage: n/oper @ h1 h2 hres with @=(+,-,*,/,@)"<<endl;
1218 return;}
1219
1220// Decode arguments
1221const char * oper = tokens[0].c_str();
1222if( oper[0]!='+' && oper[0]!='-' && oper[0]!='*' && oper[0]!='/' )
1223 {cout<<"PAWExecutor::h_oper Error: unknow operation "<<oper<<endl;
1224 return;}
1225string h1name = tokens[1];
1226string h2name = tokens[2];
1227string h3name = tokens[3];
1228
1229// Get objects
1230NamedObjMgr omg;
1231AnyDataObj* mobjh1 = omg.GetObj(h1name);
1232AnyDataObj* mobjh2 = omg.GetObj(h2name);
1233if( mobjh1==NULL || mobjh2==NULL )
1234 {cout<<"PAWExecutor::h_oper Error: unknow object(s) "<<h1name<<" or "<<h2name<<endl;
1235 return;}
1236AnyDataObj* mobjh3 = omg.GetObj(h3name);
1237
1238// Operations on HProf, only + is working
1239if( dynamic_cast<HProf*>(mobjh1) != NULL ) {
1240 if( oper[0]!='+' )
1241 { cout<<"PAWExecutor::h_oper Error: operation "<<oper
1242 <<" not implemented for HProf"<<endl; return;}
1243 if( dynamic_cast<HProf*>(mobjh2) == NULL )
1244 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h2\n"
1245 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh2).name()<<endl;
1246 return;}
1247 HProf* h1 =(HProf*) mobjh1;
1248 HProf* h2 =(HProf*) mobjh2;
1249 if( h1->NBins() != h2->NBins() )
1250 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h2 "
1251 <<h1->NBins()<<" "<<h2->NBins()<<endl; return;}
1252 HProf* h3 = NULL;
1253 if( mobjh3 == NULL ) // l'objet n'existe pas, on le cree
1254 {h3 = new HProf(*h1); h3->Zero(); omg.AddObj(h3,h3name); mobjh3 = omg.GetObj(h3name);}
1255 h3 = dynamic_cast<HProf*>(mobjh3);
1256 if(h3 == NULL) // ce n'est pas un HProf
1257 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h3\n"
1258 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh3).name()<<endl; return;}
1259 if(h1->NBins() != h3->NBins())
1260 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h3 "
1261 <<h1->NBins()<<" "<<h3->NBins()<<endl; return;}
1262 *h3 = *h1 + *h2;
1263 h3->UpdateHisto();
1264
1265// Operations on Histo
1266} else if( dynamic_cast<Histo*>(mobjh1) != NULL ) {
1267 if( dynamic_cast<Histo*>(mobjh2) == NULL )
1268 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h2\n"
1269 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh2).name()<<endl;
1270 return;}
1271 Histo* h1 =(Histo*) mobjh1;
1272 Histo* h2 =(Histo*) mobjh2;
1273 if( h1->NBins() != h2->NBins() )
1274 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h2 "
1275 <<h1->NBins()<<" "<<h2->NBins()<<endl; return;}
1276 Histo* h3 = NULL;
1277 if( mobjh3 == NULL ) // l'objet n'existe pas, on le cree
1278 {h3 = new Histo(*h1); h3->Zero(); omg.AddObj(h3,h3name); mobjh3 = omg.GetObj(h3name);}
1279 h3 = dynamic_cast<Histo*>(mobjh3);
1280 if(h3 == NULL) // ce n'est pas un Histo
1281 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h3\n"
1282 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh3).name()<<endl; return;}
1283 if(h1->NBins() != h3->NBins())
1284 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h3 "
1285 <<h1->NBins()<<" "<<h3->NBins()<<endl; return;}
1286 if( oper[0]=='+') *h3 = *h1 + *h2;
1287 else if( oper[0]=='-') *h3 = *h1 - *h2;
1288 else if( oper[0]=='*') *h3 = *h1 * *h2;
1289 else if( oper[0]=='/') *h3 = *h1 / *h2;
1290
1291// Operations on Histo2D
1292} else if( dynamic_cast<Histo2D*>(mobjh1) != NULL ) {
1293 if( dynamic_cast<Histo2D*>(mobjh2) == NULL )
1294 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h2\n"
1295 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh2).name()<<endl;
1296 return;}
1297 Histo2D* h1 =(Histo2D*) mobjh1;
1298 Histo2D* h2 =(Histo2D*) mobjh2;
1299 if( h1->NBinX() != h2->NBinX() || h1->NBinY() != h2->NBinY() )
1300 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h2 "
1301 <<h1->NBinX()<<","<<h1->NBinY()<<" "
1302 <<h2->NBinX()<<","<<h2->NBinY()<<endl; return;}
1303 Histo2D* h3 = NULL;
1304 if( mobjh3 == NULL ) // l'objet n'existe pas, on le cree
1305 {h3 = new Histo2D(*h1); h3->Zero(); omg.AddObj(h3,h3name); mobjh3 = omg.GetObj(h3name);}
1306 h3 = dynamic_cast<Histo2D*>(mobjh3);
1307 if(h3 == NULL) // ce n'est pas un Histo2D
1308 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h3\n"
1309 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh3).name()<<endl; return;}
1310 if( h1->NBinX() != h3->NBinX() || h1->NBinY() != h3->NBinY() )
1311 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h3 "
1312 <<h1->NBinX()<<","<<h1->NBinY()<<" "
1313 <<h3->NBinX()<<","<<h3->NBinY()<<endl; return;}
1314 if( oper[0]=='+') *h3 = *h1 + *h2;
1315 else if( oper[0]=='-') *h3 = *h1 - *h2;
1316 else if( oper[0]=='*') *h3 = *h1 * *h2;
1317 else if( oper[0]=='/') *h3 = *h1 / *h2;
1318
1319// Operations on Vector
1320} else if( dynamic_cast<Vector*>(mobjh1) != NULL ) {
1321 if( dynamic_cast<Vector*>(mobjh2) == NULL )
1322 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h2\n"
1323 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh2).name()<<endl;
1324 return;}
1325 Vector* h1 =(Vector*) mobjh1;
1326 Vector* h2 =(Vector*) mobjh2;
1327 if( h1->NElts() != h2->NElts() )
1328 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h2 "
1329 <<h1->NElts()<<" "<<h2->NElts()<<endl; return;}
1330 Vector* h3 = NULL;
1331 if( mobjh3 == NULL ) { // l'objet n'existe pas, on le cree
1332#ifdef SANS_EVOLPLANCK
1333 h3 = new Vector(*h1);
1334#else
1335 h3 = new Vector(*h1,false);
1336#endif
1337 *h3 = 0.; omg.AddObj(h3,h3name); mobjh3 = omg.GetObj(h3name);
1338 }
1339 h3 = dynamic_cast<Vector*>(mobjh3);
1340 if(h3 == NULL) // ce n'est pas un Vector
1341 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h3\n"
1342 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh3).name()<<endl; return;}
1343 if(h1->NElts() != h3->NElts())
1344 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h3 "
1345 <<h1->NElts()<<" "<<h3->NElts()<<endl; return;}
1346 if( oper[0]=='+') *h3 = *h1 + *h2;
1347 else if( oper[0]=='-') *h3 = *h1 - *h2;
1348#ifdef SANS_EVOLPLANCK
1349 else if(oper[0]=='*' || oper[0]=='/')
1350 cout<<"PAWExecutor::h_oper Error: operation "<<oper[0]
1351 <<" not implemented for Vector in Peida"<<endl;
1352#else
1353 else if( oper[0]=='*') {h3->Clone(*h1); h3->MulElt(*h2,*h3);}
1354 else if( oper[0]=='/') {h3->Clone(*h1); h3->DivElt(*h2,*h3,false,true);}
1355#endif
1356
1357// Operations on Matrix
1358} else if( dynamic_cast<Matrix*>(mobjh1) != NULL ) {
1359 if( dynamic_cast<Matrix*>(mobjh2) == NULL )
1360 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h2\n"
1361 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh2).name()<<endl;
1362 return;}
1363 Matrix* h1 =(Matrix*) mobjh1;
1364 Matrix* h2 =(Matrix*) mobjh2;
1365 if( h1->NRows() != h2->NRows() || h1->NCol() != h2->NCol() )
1366 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h2 "
1367 <<h1->NRows()<<","<<h1->NCol()<<" "
1368 <<h2->NRows()<<","<<h2->NCol()<<endl; return;}
1369 Matrix* h3 = NULL;
1370 if( mobjh3 == NULL ) { // l'objet n'existe pas, on le cree
1371#ifdef SANS_EVOLPLANCK
1372 h3 = new Matrix(*h1);
1373#else
1374 h3 = new Matrix(*h1,false);
1375#endif
1376 *h3 = 0.; omg.AddObj(h3,h3name); mobjh3 = omg.GetObj(h3name);
1377 }
1378 h3 = dynamic_cast<Matrix*>(mobjh3);
1379 if(h3 == NULL) // ce n'est pas un Matrix
1380 {cout<<"PAWExecutor::h_oper Error: type mismatch between h1 and h3\n"
1381 <<typeid(*mobjh1).name()<<" , "<<typeid(*mobjh3).name()<<endl; return;}
1382 if( h1->NRows() != h3->NRows() || h1->NCol() != h3->NCol() )
1383 {cout<<"PAWExecutor::h_oper Error: size mismatch between h1, h3 "
1384 <<h1->NRows()<<","<<h1->NCol()<<" "
1385 <<h3->NRows()<<","<<h3->NCol()<<endl; return;}
1386 if( oper[0]=='+') *h3 = *h1 + *h2;
1387 else if( oper[0]=='-') *h3 = *h1 - *h2;
1388#ifdef SANS_EVOLPLANCK
1389 else if(oper[0]=='*' || oper[0]=='/')
1390 cout<<"PAWExecutor::h_oper Error: operation "<<oper[0]
1391 <<" not implemented for Vector in Peida"<<endl;
1392#else
1393 else if( oper[0]=='*') {h3->Clone(*h1); h3->MulElt(*h2,*h3);}
1394 else if( oper[0]=='/') {h3->Clone(*h1); h3->DivElt(*h2,*h3,false,true);}
1395#endif
1396
1397// Doesn't work for other objects
1398} else {
1399 cout<<"PAWExecutor::h_oper Error: not implemented for "
1400 <<typeid(*mobjh1).name()<<endl;
1401 return;
1402}
1403
1404return;
1405}
1406
1407/* methode */
1408void PAWExecutor::h_plot_2d(vector<string>& tokens)
1409// plot for 2D histogramme: plot histo, bandx/y, slicex/y or projx/y
1410{
1411if(tokens.size()<1)
1412 {cout<<"Usage: h/plot/2d nameh2d to_plot [n/s] [dopt]"<<endl; return;}
1413string proj = "h"; if(tokens.size()>1) proj = tokens[1];
1414NamedObjMgr omg;
1415AnyDataObj* mobj = omg.GetObj(tokens[0]);
1416if(mobj==NULL)
1417 {cout<<"PAWExecutor::h_plot_2d Error: unknow object"<<tokens[0]<<endl;
1418 return;}
1419Histo2D* h2 = dynamic_cast<Histo2D*>(mobj);
1420if(!h2)
1421 {cout<<"PAWExecutor::h_plot_2d Error: "<<tokens[0]<<" not an Histo2D"<<endl;
1422 return;}
1423
1424Histo* h1p = NULL; string nametoplot = "/autoc/h_plot_2d_h1";
1425
1426string dopt = ""; if(tokens.size()>=3) dopt = tokens[2];
1427if(proj == "show") {
1428 h2->ShowProj();
1429 h2->ShowBand(2);
1430 h2->ShowSli(2);
1431 return;
1432} else if(proj == "h") {
1433 nametoplot = tokens[0];
1434} else if(proj == "px") {
1435 if((h1p=h2->HProjX())) {Histo* h1=new Histo(*h1p); omg.AddObj(h1,nametoplot);}
1436 else {h2->ShowProj(); return;}
1437} else if(proj == "py") {
1438 if((h1p=h2->HProjY())) {Histo* h1=new Histo(*h1p); omg.AddObj(h1,nametoplot);}
1439 else {h2->ShowProj(); return;}
1440} else {
1441 if(tokens.size()<3)
1442 {cout<<"Usage: h/plot/2d nameh2d bx/by/sx/sy n [dopt]"<<endl; return;}
1443 int_4 n = atoi(tokens[2].c_str());
1444 dopt = ""; if(tokens.size()>=4) dopt = tokens[3];
1445 if(proj == "bx") {
1446 if((h1p=h2->HBandX(n))) {Histo* h1=new Histo(*h1p); omg.AddObj(h1,nametoplot);}
1447 else {h2->ShowBand(); return;}
1448 } else if(proj == "by") {
1449 if((h1p=h2->HBandY(n))) {Histo* h1=new Histo(*h1p); omg.AddObj(h1,nametoplot);}
1450 else {h2->ShowBand(); return;}
1451 } else if(proj == "sx") {
1452 if((h1p=h2->HSliX(n))) {Histo* h1=new Histo(*h1p); omg.AddObj(h1,nametoplot);}
1453 else {h2->ShowSli(); return;}
1454 } else if(proj == "sy") {
1455 if((h1p=h2->HSliY(n))) {Histo* h1=new Histo(*h1p); omg.AddObj(h1,nametoplot);}
1456 else {h2->ShowSli(); return;}
1457 }
1458}
1459
1460omg.DisplayObj(nametoplot,dopt);
1461}
1462
1463/* methode */
1464int_4 PAWExecutor::decodepawstring(string tokens,string& nameobj
1465 ,string& xexp,string& yexp,string& zexp)
1466// Decodage general de "nameobj.xexp"
1467// "nameobj.yexp%xexp"
1468// "nameobj.zexp%yexp%xexp"
1469// Return: nombre de variables trouvees, -1 si probleme
1470{
1471nameobj = ""; xexp= ""; yexp= ""; zexp= "";
1472
1473int_4 lt = (int) tokens.length();
1474if(lt<=0) return -1;
1475
1476// decodage de la chaine de type PAW.
1477char *str = new char[lt+2];
1478strcpy(str,tokens.c_str()); strip(str,'B',' '); lt = strlen(str);
1479//cout<<"chaine1["<<lt<<"] :"<<str<<":"<<endl;
1480char *c[3] = {NULL,NULL,NULL};
1481int_4 i, np=0; bool namefound = false;
1482for(i=0;i<lt;i++) {
1483 if(!namefound && str[i]=='.') {
1484 str[i]='\0';
1485 namefound=true;
1486 c[np] = str+i+1; np++;
1487 }
1488 if( namefound && str[i]=='%') {
1489 str[i]='\0';
1490 if(np<3) {c[np] = str+i+1; np++;}
1491 }
1492}
1493//cout<<"chaine2 :"; for(i=0;i<lt;i++) cout<<str[i]; cout<<":"<<endl;
1494
1495// Remplissage du nom et des variables
1496nameobj = str;
1497if(np==1) xexp=c[0];
1498if(np==2) {yexp=c[0]; xexp=c[1];}
1499if(np==3) {zexp=c[0]; yexp=c[1]; xexp=c[2];}
1500//cout<<"pawstring str,c[0-2] "<<str<<" "<<c[0]<<" "<<c[1]<<" "<<c[2]<<endl;
1501delete [] str;
1502
1503// Comptage des variables
1504np = -1;
1505if(nameobj.length()>0)
1506 {np = 0; if(xexp.length()>0)
1507 {np++; if(yexp.length()>0)
1508 {np++; if(zexp.length()>0) np++;}}}
1509//cout<<"pawstring["<<np<<"] name="<<nameobj
1510// <<" xexp="<<xexp<<" yexp="<<yexp<<" zexp="<<zexp<<endl;
1511return np;
1512}
1513
1514/* methode */
1515void PAWExecutor::h_put_vec(vector<string>& tokens)
1516// Pour remplir un histo avec le contenu d'un vecteur ou d'une matrice
1517// L'histogramme doit deja exister. Le nombre de bins remplit
1518// depend des tailles respectives des 2 objets.
1519// tokens[2] = "cont" : on remplit les valeurs de l'histogramme (ou "!")
1520// = "err2" : on remplit les erreurs de l'histogramme
1521{
1522if(tokens.size()<2)
1523 {cout<<"Usage: h/put_vec namehisto namevector"<<endl;
1524 return;}
1525string toput = "cont"; if(tokens.size()>2) toput = tokens[2];
1526if(toput=="!") toput = "cont";
1527if(toput!="cont" && toput!="err2")
1528 {cout<<"PAWExecutor::h_put_vec Error: unknow filling "<<toput<<endl;
1529 return;}
1530string hname = tokens[0];
1531string vname = tokens[1];
1532
1533// Get objects
1534NamedObjMgr omg;
1535AnyDataObj* mobjh = omg.GetObj(hname);
1536AnyDataObj* mobjv = omg.GetObj(vname);
1537if( mobjh==NULL || mobjv==NULL )
1538 {cout<<"PAWExecutor::h_put_vec Error: unknow object(s) "<<hname<<" or "<<vname<<endl;
1539 return;}
1540
1541// Fill Histo from Vector
1542if(typeid(*mobjh) == typeid(Histo) && typeid(*mobjv) == typeid(Vector)) {
1543 Histo* h = dynamic_cast<Histo*>(mobjh);
1544 Vector* v = dynamic_cast<Vector*>(mobjv);
1545 if(toput=="cont") h->PutValue(*v);
1546 else if(toput=="err2") h->PutError2(*v);
1547
1548// Fill Histo2D from Matrix
1549} else if(typeid(*mobjh) == typeid(Histo2D) && typeid(*mobjv) == typeid(Matrix)) {
1550 Histo2D* h = dynamic_cast<Histo2D*>(mobjh);
1551 Matrix* v = dynamic_cast<Matrix*>(mobjv);
1552 if(toput=="cont") h->PutValue(*v);
1553 else if(toput=="err2") h->PutError2(*v);
1554
1555} else {
1556 cout<<"PAWExecutor::h_put_vec Error: type mismatch between histogram and vector/matrix\n"
1557 <<typeid(*mobjh).name()<<" , "<<typeid(*mobjv).name()<<endl;
1558 return;
1559}
1560
1561}
1562
1563/* methode */
1564void PAWExecutor::h_get_vec(vector<string>& tokens)
1565// Pour copier un histo dans un vecteur ou une matrice
1566// Si le vecteur (matrice) n'existe pas, il est cree.
1567// Le vecteur ou la matrice est re-dimensionne correctement.
1568// tokens[2] = "cont" : on copie les valeurs de l'histogramme (ou "!")
1569// = "err2" : on copie les erreurs de l'histogramme
1570// = "absc" : on copie les erreurs de l'histogramme (1D only)
1571// tokens[3[,4]] = show : show available projections for Histo2D
1572// px : get X projection
1573// py : get Y projection
1574// bx n : get X band number n
1575// by n : get Y band number n
1576// sx n : get X slice number n
1577// sy n : get Y slice number n
1578{
1579if(tokens.size()<2)
1580 {cout<<"Usage: h/get_vec namehisto namevector"<<endl;
1581 return;}
1582string toget = "cont"; if(tokens.size()>2) toget = tokens[2];
1583if(toget=="!") toget = "cont";
1584if(toget!="cont" && toget!="err2" && toget!="absc")
1585 {cout<<"PAWExecutor::h_get_vec Error: unknow filling "<<toget<<endl;
1586 return;}
1587string proj = "h"; if(tokens.size()>3) proj = tokens[3];
1588int_4 nproj = -1; if(tokens.size()>4) nproj = atoi(tokens[4].c_str());
1589string hname = tokens[0];
1590string vname = tokens[1];
1591
1592// Get objects
1593NamedObjMgr omg;
1594AnyDataObj* mobjh = omg.GetObj(hname);
1595AnyDataObj* mobjv = omg.GetObj(vname);
1596if( mobjh==NULL)
1597 {cout<<"PAWExecutor::h_put_vec Error: unknow object(s) "<<hname<<endl;
1598 return;}
1599
1600// Copy Histo/Histo2D/Projection to Vector/Matrix
1601Histo* h = dynamic_cast<Histo*>(mobjh);
1602Histo2D* h2 = dynamic_cast<Histo2D*>(mobjh);
1603if( h != NULL ) { // Histo ou HProf
1604 h->UpdateHisto(); // pour le cas ou c'est un HProf
1605 Vector* v = NULL;
1606 if(mobjv==NULL) // le vecteur n'existe pas
1607 {v = new Vector(1); omg.AddObj(v,vname); mobjv = omg.GetObj(vname);}
1608 if(typeid(*mobjv) != typeid(Vector))
1609 {cout<<"PAWExecutor::h_get_vec Error: type mismatch between Histo/HProf and vector\n"
1610 <<typeid(*mobjv).name()<<endl; return;}
1611 v = dynamic_cast<Vector*>(mobjv);
1612 if(toget=="cont") h->GetValue(*v);
1613 else if(toget=="err2") h->GetError2(*v);
1614 else if(toget=="absc") h->GetAbsc(*v);
1615
1616} else if( h2 != NULL) { // Histo2D
1617
1618 if(proj == "h") { // On veut les valeurs de l'histogramme 2D
1619 Matrix* v = NULL;
1620 if(mobjv==NULL) // la matrice n'existe pas
1621 {v = new Matrix(1,1); omg.AddObj(v,vname); mobjv = omg.GetObj(vname);}
1622 if(typeid(*mobjv) != typeid(Matrix))
1623 {cout<<"PAWExecutor::h_get_vec Error: type mismatch between Histo2D and matrix\n"
1624 <<typeid(*mobjv).name()<<endl; return;}
1625 v = dynamic_cast<Matrix*>(mobjv);
1626 if(toget=="cont") h2->GetValue(*v);
1627 else if(toget=="err2") h2->GetError2(*v);
1628 else
1629 {cout<<"PAWExecutor::h_get_vec Error: option "<<toget<<" not valid for Histo2D"<<endl;
1630 return;}
1631
1632 } else { // On veut les valeurs d'une projection de l'histo 2D
1633 h = NULL;
1634 if(proj == "show") {h2->ShowProj(); h2->ShowBand(2); h2->ShowSli(2);}
1635 else if(proj == "px") h = h2->HProjX();
1636 else if(proj == "py") h = h2->HProjY();
1637 else if(proj == "bx") h = h2->HBandX(nproj);
1638 else if(proj == "by") h = h2->HBandY(nproj);
1639 else if(proj == "sx") h = h2->HSliX(nproj);
1640 else if(proj == "sy") h = h2->HSliY(nproj);
1641 if(h==NULL)
1642 {cout<<"PAWExecutor::h_get_vec Error: unknown projection "<<proj
1643 <<" number "<<nproj<<endl; return;}
1644 Vector* v = NULL;
1645 if(mobjv==NULL) // le vecteur n'existe pas
1646 {v = new Vector(1); omg.AddObj(v,vname); mobjv = omg.GetObj(vname);}
1647 if(typeid(*mobjv) != typeid(Vector))
1648 {cout<<"PAWExecutor::h_get_vec Error: type mismatch between Histo2D "<<proj
1649 <<" and vector\n"<<typeid(*mobjv).name()<<endl; return;}
1650 v = dynamic_cast<Vector*>(mobjv);
1651 if(toget=="cont") h->GetValue(*v);
1652 else if(toget=="err2") h->GetError2(*v);
1653 else if(toget=="absc") h->GetAbsc(*v);
1654 }
1655
1656} else {
1657 cout<<"PAWExecutor::h_get_vec Error: type mismatch for histogram\n"
1658 <<typeid(*mobjh).name()<<endl;
1659 return;
1660}
1661
1662}
1663
1664/* methode */
1665void PAWExecutor::h_copy(vector<string>& tokens)
1666// Pour copier un object dans un object
1667// objects are Vector,Matrix,Histo,Histo2D
1668// h/copy namefrom nametocopy [i1[:i2]] [j1[:j2]] [ic1[:jc1]]
1669// copy obj1Dfrom(i1->i2) into obj1Dto(ic1->)
1670// copy obj2Dfrom(i1,j1->i2,j2) into obj2Dto(ic1,jc1->)
1671// Attention: elements depuis i1 jusqu'a i2 COMPRIS
1672// Si obj1Dto n'existe pas, il est cree avec une taille i2-i1+1
1673// ou obj2Dto avec une taille i2-i1+1,j2-j1+1
1674// Le contenu de obj?Dfrom adresse est surecrit
1675// Le contenu de obj?Dfrom non adresse reste le meme
1676{
1677int_4 tks = tokens.size();
1678if(tks<2)
1679 {cout<<"Usage: h_copy namefrom nametocopy [i1[:i2]] [j1[:j2]] [ic1[:jc1]]"<<endl;
1680 return;}
1681
1682NamedObjMgr omg;
1683AnyDataObj* mobjv1 = omg.GetObj(tokens[0]);
1684if( mobjv1==NULL)
1685 {cout<<"PAWExecutor::h_copy Error: unknow object "<<tokens[0]<<endl;
1686 return;}
1687AnyDataObj* mobjv2 = omg.GetObj(tokens[1]);
1688
1689int_4 i1=0,i2=0, j1=0,j2=0, ic1=0,jc1=0, i,ii, j,jj, nx1,ny1, nx2,ny2;
1690
1691// Cas d'un Vector
1692if(typeid(*mobjv1) == typeid(Vector)) {
1693 Vector* v1 = dynamic_cast<Vector*>(mobjv1); nx1 = (int_4)v1->NElts();
1694 i2=nx1-1;
1695 if(tks>2) if(tokens[2]!="!") sscanf(tokens[2].c_str(),"%d:%d",&i1,&i2);
1696 if(i1<0) i1=0; if(i2>=nx1) i2=nx1-1;
1697 if(i2<i1)
1698 {cout<<"PAWExecutor::h_copy Error: wrong range ["<<i1<<":"<<i2<<"] for NElts="
1699 <<nx1<<endl; return;}
1700 Vector* v2 = NULL;
1701 if(mobjv2==NULL)
1702 {v2 = new Vector(i2-i1+1); omg.AddObj(v2,tokens[1]); mobjv2 = omg.GetObj(tokens[1]);}
1703 if(typeid(*mobjv2) != typeid(Vector))
1704 {cout<<"PAWExecutor::h_copy Error: type mismatch for Vector "
1705 <<typeid(*mobjv2).name()<<endl; return;}
1706 v2 = dynamic_cast<Vector*>(mobjv2); nx2 = (int_4)v2->NElts();
1707 if(tks>4) if(tokens[4]!="!") sscanf(tokens[4].c_str(),"%d",&ic1);
1708 if(ic1<0) ic1=0;
1709 if(ic1>=nx2)
1710 {cout<<"PAWExecutor::h_copy Error: wrong pointer to copy ["<<ic1<<"] for NElts="
1711 <<nx2<<endl; return;}
1712 cout<<"copy "<<tokens[0]<<"("<<i1<<":"<<i2<<") to "<<tokens[1]<<"("<<ic1<<"->...)"<<endl;
1713 for(i=i1,ii=ic1; i<=i2 && i<nx1 && ii<nx2; i++,ii++) (*v2)(ii) = (*v1)(i);
1714 return;
1715}
1716
1717// Cas d'un Histo
1718if(typeid(*mobjv1) == typeid(Histo)) {
1719 Histo* v1 = dynamic_cast<Histo*>(mobjv1); nx1 = (int_4)v1->NBins();
1720 i2=nx1-1;
1721 if(tks>2) if(tokens[2]!="!") sscanf(tokens[2].c_str(),"%d:%d",&i1,&i2);
1722 if(i1<0) i1=0; if(i2>=nx1) i2=nx1-1;
1723 if(i2<i1)
1724 {cout<<"PAWExecutor::h_copy Error: wrong range ["<<i1<<":"<<i2<<"] for NBins="
1725 <<nx1<<endl; return;}
1726 Histo* v2 = NULL;
1727 if(mobjv2==NULL)
1728 {v2 = new Histo(0.,(r_8)(i2-i1+1),i2-i1+1);
1729 omg.AddObj(v2,tokens[1]); mobjv2 = omg.GetObj(tokens[1]);}
1730 if(typeid(*mobjv2) != typeid(Histo))
1731 {cout<<"PAWExecutor::h_copy Error: type mismatch for Histo "
1732 <<typeid(*mobjv2).name()<<endl; return;}
1733 v2 = dynamic_cast<Histo*>(mobjv2); nx2 = (int_4)v2->NBins();
1734 if(v1->HasErrors() && !(v2->HasErrors())) v2->Errors();
1735 if(tks>4) if(tokens[4]!="!") sscanf(tokens[4].c_str(),"%d",&ic1);
1736 if(ic1<0) ic1=0;
1737 if(ic1>=nx2)
1738 {cout<<"PAWExecutor::h_copy Error: wrong pointer to copy ["<<ic1<<"] for NBins="
1739 <<nx2<<endl; return;}
1740 cout<<"copy "<<tokens[0]<<"("<<i1<<":"<<i2<<") to "<<tokens[1]<<"("<<ic1<<"->...)"<<endl;
1741 for(i=i1,ii=ic1; i<=i2 && i<nx1 && ii<nx2; i++,ii++)
1742 {(*v2)(ii) = (*v1)(i); if(v1->HasErrors()) v2->Error2(ii) = v1->Error2(i);}
1743 return;
1744}
1745
1746// Cas d'une Matrix
1747if(typeid(*mobjv1) == typeid(Matrix)) {
1748 Matrix* v1 = dynamic_cast<Matrix*>(mobjv1); nx1=v1->NRows(); ny1=v1->NCol();
1749 i2=nx1-1; j2=ny1-1;
1750 if(tks>2) if(tokens[2]!="!") sscanf(tokens[2].c_str(),"%d:%d",&i1,&i2);
1751 if(tks>3) if(tokens[3]!="!") sscanf(tokens[3].c_str(),"%d:%d",&j1,&j2);
1752 if(i1<0) i1=0; if(i2>=nx1) i2=nx1-1; if(j1<0) j1=0; if(j2>=ny1) j2=ny1-1;
1753 if(i2<i1 || j2<j1)
1754 {cout<<"PAWExecutor::h_copy Error: wrong range ["<<i1<<":"<<i2
1755 <<"] , ["<<j1<<":"<<j2<<"] for NRows,NCol="
1756 <<nx1<<" , "<<ny1<<endl; return;}
1757 Matrix* v2 = NULL;
1758 if(mobjv2==NULL)
1759 {v2 = new Matrix(i2-i1+1,j2-j1+1);
1760 omg.AddObj(v2,tokens[1]); mobjv2 = omg.GetObj(tokens[1]);}
1761 if(typeid(*mobjv2) != typeid(Matrix))
1762 {cout<<"PAWExecutor::h_copy Error: type mismatch for Matrix "
1763 <<typeid(*mobjv2).name()<<endl; return;}
1764 v2 = dynamic_cast<Matrix*>(mobjv2); nx2=v2->NRows(); ny2=v2->NCol();
1765 if(tks>4) if(tokens[4]!="!") sscanf(tokens[4].c_str(),"%d:%d",&ic1,&jc1);
1766 if(ic1<0) ic1=0; if(jc1<0) jc1=0;
1767 if(ic1>=nx2 || jc1>=ny2)
1768 {cout<<"PAWExecutor::h_copy Error: wrong pointer to copy ["<<ic1<<","<<jc1
1769 <<"] for NRows,NCol="<<nx2<<","<<ny2<<endl; return;}
1770 cout<<"copy "<<tokens[0]<<"("<<i1<<":"<<i2<<","<<j1<<":"<<j2
1771 <<") to "<<tokens[1]<<"("<<ic1<<","<<jc1<<"->...)"<<endl;
1772 for(i=i1,ii=ic1; i<=i2 && i<nx1 && ii<nx2; i++,ii++)
1773 for(j=j1,jj=jc1; j<=j2 && j<ny1 && jj<ny2; j++,jj++) (*v2)(ii,jj) = (*v1)(i,j);
1774 return;
1775}
1776
1777// Cas d'un Histo2D
1778if(typeid(*mobjv1) == typeid(Histo2D)) {
1779 Histo2D* v1 = dynamic_cast<Histo2D*>(mobjv1); nx1=v1->NBinX(); ny1=v1->NBinY();
1780 i2=nx1-1; j2=ny1-1;
1781 if(tks>2) if(tokens[2]!="!") sscanf(tokens[2].c_str(),"%d:%d",&i1,&i2);
1782 if(tks>3) if(tokens[3]!="!") sscanf(tokens[3].c_str(),"%d:%d",&j1,&j2);
1783 if(i1<0) i1=0; if(i2>=nx1) i2=nx1-1; if(j1<0) j1=0; if(j2>=ny1) j2=ny1-1;
1784 if(i2<i1 || j2<j1)
1785 {cout<<"PAWExecutor::h_copy Error: wrong range ["<<i1<<":"<<i2
1786 <<"] , ["<<j1<<":"<<j2<<"] for NBinX,NBinY="
1787 <<nx1<<" , "<<ny1<<endl; return;}
1788 Histo2D* v2 = NULL;
1789 if(mobjv2==NULL)
1790 {v2 = new Histo2D(0.,(r_8)(i2-i1+1),i2-i1+1,0.,(r_8)(j2-j1+1),j2-j1+1);
1791 omg.AddObj(v2,tokens[1]); mobjv2 = omg.GetObj(tokens[1]);}
1792 if(typeid(*mobjv2) != typeid(Histo2D))
1793 {cout<<"PAWExecutor::h_copy Error: type mismatch for Histo2D"
1794 <<typeid(*mobjv2).name()<<endl; return;}
1795 v2 = dynamic_cast<Histo2D*>(mobjv2); nx2=v2->NBinX(); ny2=v2->NBinY();
1796 if(v1->HasErrors() && !(v2->HasErrors())) v2->Errors();
1797 if(tks>4) sscanf(tokens[4].c_str(),"%d:%d",&ic1,&jc1);
1798 if(ic1<0) ic1=0; if(jc1<0) jc1=0;
1799 if(ic1>=nx2 || jc1>=ny2)
1800 {cout<<"PAWExecutor::h_copy Error: wrong pointer to copy ["<<ic1<<","<<jc1
1801 <<"] for NBinX,NBinY="<<nx2<<","<<ny2<<endl; return;}
1802 cout<<"copy "<<tokens[0]<<"("<<i1<<":"<<i2<<","<<j1<<":"<<j2
1803 <<") to "<<tokens[1]<<"("<<ic1<<","<<jc1<<"->...)"<<endl;
1804 for(i=i1,ii=ic1; i<=i2 && i<nx1 && ii<nx2; i++,ii++)
1805 for(j=j1,jj=jc1; j<=j2 && j<ny1 && jj<ny2; j++,jj++)
1806 {(*v2)(ii,jj) = (*v1)(i,j); if(v1->HasErrors()) v2->Error2(ii,jj) = v1->Error2(i,j);}
1807 return;
1808}
1809
1810// Cas non prevu
1811cout<<"PAWExecutor::h_copy Error: type mismatch for Vector/Matrix/Histo/Histo2D\n"
1812 <<typeid(*mobjv1).name()<<endl;
1813return;
1814}
1815
1816/* methode */
1817void PAWExecutor::h_set(string dum,vector<string>& tokens)
1818// Mise de valeurs/erreurs d'un histo 1D ou 2D a une valeur donnee
1819// dum = err : on change les valeurs des erreurs
1820// cont : on change les valeurs du contenu des bins
1821// tokens[0] : nom de l'histogramme
1822// tokens[1] : valeur de remplacement
1823// tokens[2-3] : i1:i2 j1:j2 intervalle des bins qui doivent etre remplace
1824// : v v1:v2 remplacement des bins ayant une valeur dans cet intervalle
1825// : e e1:e2 remplacement des bins ayant une erreur dans cet intervalle
1826{
1827int_4 tks = tokens.size();
1828if(tks<3)
1829 {cout<<"Usage: h/set/[err,cont] namehisto setvalue i1[:i2] [j1[:j2]]\n"
1830 <<" v v1:v2\n"
1831 <<" e e1:e2"
1832 <<endl; return;}
1833
1834// Decodage arguments
1835bool replerr = false; if(dum=="err") replerr = true;
1836r_8 setval = atof(tokens[1].c_str());
1837int_4 testcont=0; if(tokens[2]=="v") testcont=1; if(tokens[2]=="e") testcont=2;
1838int_4 i1=-1, i2=-1, j1=-1, j2=-1;
1839r_8 v1=0., v2=0.;
1840if(testcont==0) {
1841 sscanf(tokens[2].c_str(),"%d:%d",&i1,&i2);
1842 if(tks>3) sscanf(tokens[3].c_str(),"%d:%d",&j1,&j2);
1843} else {
1844 if(tks<=3)
1845 {cout<<"PAWExecutor::h_set Error: h/set/... v v1:v2, please give v1:v2"<<endl;
1846 return;}
1847 sscanf(tokens[3].c_str(),"%lf:%lf",&v1,&v2);
1848}
1849
1850if(replerr) {if(setval<0.) setval = 0.; else setval *= setval;}
1851if(testcont!=0 && v2<v1)
1852 {cout<<"PAWExecutor::h_set Error: bad value range="<<v1<<","<<v2<<endl; return;}
1853
1854// identification objet
1855NamedObjMgr omg;
1856AnyDataObj* mobj = omg.GetObj(tokens[0]);
1857if( mobj==NULL)
1858 {cout<<"PAWExecutor::h_set Error: unknow object "<<tokens[0]<<endl;
1859 return;}
1860
1861// Traitement
1862if(typeid(*mobj) == typeid(Histo)) { // Histo 1D
1863 Histo* h = dynamic_cast<Histo*>(mobj);
1864 if( (replerr || testcont==2) && !(h->HasErrors()) )
1865 {cout<<"PAWExecutor::h_set Error: histogram has no errors"<<endl; return;}
1866 if(testcont!=0) {i1=0; i2=h->NBins()-1;}
1867 if(i1<0 || i1>=h->NBins())
1868 {cout<<"PAWExecutor::h_set Error: bad bin range i1="<<i1<<endl; return;}
1869 if(i2<i1) i2=i1; if(i2>=h->NBins()) i2=h->NBins()-1;
1870 for(int_4 i=i1;i<=i2;i++) {
1871 bool change = true;
1872 if(testcont==1) {if((*h)(i)<v1 || (*h)(i)>v2) change = false;}
1873 else if(testcont==2) {if(h->Error(i)<v1 || h->Error(i)>v2) change = false;}
1874 if(!change) continue;
1875 if(replerr) h->Error2(i) = setval; else (*h)(i) = setval;
1876 }
1877
1878} else if(typeid(*mobj) == typeid(Histo2D)) { // Histo 2D
1879 Histo2D* h2 = dynamic_cast<Histo2D*>(mobj);
1880 if( (replerr || testcont==2) && !(h2->HasErrors()) )
1881 {cout<<"PAWExecutor::h_set Error: histogram has no errors"<<endl; return;}
1882 if(testcont!=0) {i1=0; i2=h2->NBinX()-1;j1=0; j2=h2->NBinY()-1;}
1883 if(i1<0 || i1>=h2->NBinX() || j1<0 || j1>=h2->NBinY())
1884 {cout<<"PAWExecutor::h_set Error: bad bin range i1,j1="<<i1<<","<<j1<<endl; return;}
1885 if(i2<i1) i2=i1; if(i2>=h2->NBinX()) i2=h2->NBinX()-1;
1886 if(j2<j1) j2=j1; if(j2>=h2->NBinY()) j2=h2->NBinY()-1;
1887 for(int_4 i=i1;i<=i2;i++) for(int_4 j=j1;j<=j2;j++) {
1888 bool change = true;
1889 if(testcont==1) {if((*h2)(i,j)<v1 || (*h2)(i,j)>v2) change = false;}
1890 else if(testcont==2) {if(h2->Error(i,j)<v1 || h2->Error(i,j)>v2) change = false;}
1891 if(!change) continue;
1892 if(replerr) h2->Error2(i,j) = setval; else (*h2)(i,j) = setval;
1893 }
1894
1895} else {
1896 cout<<"PAWExecutor::h_set Error: type mismatch for Histo/Histo2D\n"
1897 <<typeid(*mobj).name()<<endl;
1898 return;
1899}
1900
1901}
1902
1903/* methode */
1904void PAWExecutor::h_err(vector<string>& tokens)
1905// Mise des erreurs d'un histo 1D ou 2D a une valeur
1906// donnee par une fonction de la valeur du bin
1907// tokens[0] : nom de l'histogramme
1908// tokens[1] : expression de la fonction
1909{
1910if(tokens.size()<2)
1911 {cout<<"Usage: h/err namehisto expr_func"<<endl; return;}
1912
1913// identification objet
1914NamedObjMgr omg;
1915AnyDataObj* mobj = omg.GetObj(tokens[0]);
1916if( mobj==NULL)
1917 {cout<<"PAWExecutor::h_err Error: unknow object "<<tokens[0]<<endl;
1918 return;}
1919
1920// Fonction
1921FILE *fip = NULL;
1922string expfunc = ""; {for(int_4 i=1;i<(int)tokens.size();i++) expfunc += tokens[i];}
1923string fname = "func1or2_pia_dl.c";
1924string func = "func1or2_pia_dl_func";
1925if((fip = fopen(fname.c_str(), "w")) == NULL)
1926 {cout<<"PAWExecutor::h_err Error: open function file "<<fname<<endl;
1927 return;}
1928fprintf(fip,"#include <math.h>\n");
1929fprintf(fip,"double %s(double x) \n{\n",func.c_str());
1930fprintf(fip,"return(%s); \n}\n", expfunc.c_str());
1931fclose(fip);
1932DlFunctionOfX f = (DlFunctionOfX) omg.GetServiceObj()->LinkFunctionFromFile(fname,func);
1933if(!f)
1934 {cout<<"PAWExecutor::h_err Error: bad function "<<func<<","<<fname<<endl;
1935 return;}
1936
1937// Traitement
1938if(typeid(*mobj) == typeid(Histo)) { // Histo 1D
1939 Histo* h = dynamic_cast<Histo*>(mobj);
1940 if(!(h->HasErrors())) h->Errors();
1941 for(int_4 i=0;i<h->NBins();i++) {
1942 r_8 x = (*h)(i);
1943 r_8 e = f(x);
1944 h->SetErr2(i,e*e);
1945 }
1946
1947} else if(typeid(*mobj) == typeid(Histo2D)) { // Histo 2D
1948 Histo2D* h2 = dynamic_cast<Histo2D*>(mobj);
1949 if(!(h2->HasErrors())) h2->Errors();
1950 for(int_4 i=0;i<h2->NBinX();i++) for(int_4 j=0;j<h2->NBinY();j++) {
1951 r_8 x = (*h2)(i,j);
1952 r_8 e = f(x);
1953 h2->Error2(i,j) = e*e;
1954 }
1955
1956} else {
1957 cout<<"PAWExecutor::h_err Error: type mismatch for Histo/Histo2D\n"
1958 <<typeid(*mobj).name()<<endl;
1959 return;
1960}
1961
1962}
Note: See TracBrowser for help on using the repository browser.