Sigbuilder block was opening two figures when asked to launch graphic window
[scilab.git] / scilab / modules / scicos_blocks / macros / Sources / CURVE_c.sci
1 //  Scicos
2 //
3 //  Copyright (C) INRIA - Masoud Najafi <masoud.najafi@inria.fr>
4 //                        Serge Steer <serge.steer@inria.fr>       1993
5 //                        Habib Jreij                              1993
6 //
7 // This program is free software; you can redistribute it and/or modify
8 // it under the terms of the GNU General Public License as published by
9 // the Free Software Foundation; either version 2 of the License, or
10 // (at your option) any later version.
11 //
12 // This program is distributed in the hope that it will be useful,
13 // but WITHOUT ANY WARRANTY; without even the implied warranty of
14 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 // GNU General Public License for more details.
16 //
17 // You should have received a copy of the GNU General Public License
18 // along with this program; if not, write to the Free Software
19 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
20 //
21 // See the file ../license.txt
22 //
23
24 function [x,y,typ]=CURVE_c(job,arg1,arg2)
25 //** 07/01/2008 : Adapted fot Scilab 5.0 by Simone Mannori 
26 x=[]; y=[]; typ=[];
27
28 select job
29
30  case "plot" then
31    standard_draw(arg1);
32
33  case "getinputs" then
34    [x,y,typ]=standard_inputs(arg1);
35
36  case "getoutputs" then
37    [x,y,typ]=standard_outputs(arg1);
38
39  case "getorigin" then
40    [x,y]=standard_origin(arg1);
41
42  case "set" then
43
44   x = arg1; 
45   model = arg1.model;
46   graphics = arg1.graphics; 
47   exprs = graphics.exprs;
48   ok = %f;
49   SaveExit = %f;
50
51   while %t do
52     Ask_again = %f; 
53     [ok,Method,xx,yy,PeriodicOption,graf,exprs] = scicos_getvalue('Spline data',['Spline"+...
54                                                   " Method (0..7)';'x';'y';'Periodic signal(y/n)?';'Launch"+...
55                                                   " graphic window(y/n)?'],list('vec',1,'vec',-1, ...
56                                                                                 'vec',-1,'str',1,'str',1),exprs)
57     if  ~ok then break;end
58     
59     if PeriodicOption=='y' | PeriodicOption=='Y' then
60       PO=1;
61     else
62       exprs(4)='n';
63       PO=0;
64     end
65     
66     mtd=int(Method);
67     if mtd<0 then
68       mtd=0
69     end;
70
71     if mtd>7 then
72       mtd=7;
73     end    
74     
75     METHOD = getmethod(mtd);
76
77     if ~Ask_again then 
78       xx=xx(:); yy=yy(:);
79       [nx,mx] = size(xx); [ny,my]=size(yy);
80       if ~((nx==ny)&(mx==my)) then
81          messagebox("Incompatible size of [x] and [y]","modal","error");
82          Ask_again = %t;
83       end
84     end
85     
86     if ~Ask_again then //+++++++++++++++++++++++++++++++++++++++
87       xy = [xx,yy];
88       [xy] = cleandata(xy); // just for sorting to be able to compare data before and after poke_point(.)
89       N= size(xy,'r');
90       exprs(5)='n';// exprs.graf='n'
91       
92       if graf=='y' | graf=='Y' then //_______Graphic editor___________
93         ipar=[N;mtd;PO];
94         rpar=[];
95
96         if (winsid() == []) then
97           curwin = 0;
98         else
99           curwin = max(winsid())+1; //** prepare a brand new win
100         end
101         //** see below in this file; "poke_point" is very similar to "edit_curv"
102         [orpar,oipar,ok] = poke_point(xy,ipar,rpar); //** HERE WE ARE +++++++++++++++++++++++++++++++++++
103         if ~ok then break;end;//  exit without save
104
105         // verifying the data change
106         N2=oipar(1);xy2=[orpar(1:N2),orpar(N2+1:2*N2)];
107         New_methhod=oipar(2);
108         DChange=%f;     
109         METHOD=getmethod(New_methhod);
110         if or(xy(:,1)<>xy2(:,1)) then, DChange=%t;end
111         if or(xy(1:N-1,2)<>xy2(1:N2-1,2)) then, DChange=%t;end
112         if (xy(N,2)<>xy2(N2,2) & (METHOD<>'periodic')) then, DChange=%t;end
113         if DChange then 
114           exprs(2)=strcat(sci2exp(xy2(:,1)))
115           exprs(3)=strcat(sci2exp(xy2(:,2)))
116         end
117
118         exprs(1)=sci2exp(New_methhod);
119         if oipar(3)==1 then,perop='y';else,perop='n';end
120         exprs(4)=perop;
121         SaveExit=%t
122        else//_____________________No graphics__________________________
123         [Xdummy,Ydummy,orpar]=Do_Spline(N,mtd,xy(:,1),xy(:,2));
124         if (METHOD=='periodic') then // periodic spline
125           xy(N,2)=xy(1,2);
126         end     
127         if (METHOD=='order 2' | METHOD=='not_a_knot'|METHOD=='periodic' | METHOD=='monotone'| METHOD=='fast' | METHOD=='clamped') then 
128           orpar=[xy(:,1);xy(:,2);orpar];                
129         else
130           if (METHOD=='zero order'|METHOD=='linear')
131             orpar=[xy(:,1);xy(:,2);]
132           end   
133         end
134         exprs(1)=sci2exp(mtd);// pour le cas methode>7 | method<0
135         oipar=[N;mtd;PO]        
136         SaveExit=%t
137       end //___________________________________________________________
138     end //++++++++++++++++++++++++++++++++++++++++++++++++++++++
139     
140     if (SaveExit) then            
141       xp=find(orpar(1:oipar(1))>=0);
142       if (xp<>[]) then 
143         model.firing=orpar(xp(1)); //first positive event
144       else  
145         model.firing=-1;
146       end
147       model.rpar=orpar
148       model.ipar=oipar
149       graphics.exprs=exprs;
150       x.model=model      
151       x.graphics=graphics
152       break
153     end
154   end
155  case 'define' then  
156   model=scicos_model()
157   xx=[0, 1, 2];yy=[10, 20, -30];  N=3;  Method=3;  PeriodicOption='y';  Graf='n'
158   model.sim=list('curve_c',4)
159   model.in=[]
160   model.out=1
161   model.rpar=[xx(:);yy(:)]
162   model.ipar=[N;Method;1]
163   model.blocktype='c'
164   model.dep_ut=[%f %t]
165   model.evtin=1
166   model.evtout=1
167   model.firing=0
168   exprs=[sci2exp(Method);sci2exp(xx);sci2exp(yy);PeriodicOption;Graf]
169   
170   gr_i=['rpar=arg1.model.rpar;n=model.ipar(1);order=model.ipar(2);';
171         'xx=rpar(1:n);yy=rpar(n+1:2*n);';
172         '[XX,YY,rpardummy]=Do_Spline(n,order,xx,yy)';
173         'xmx=max(XX);xmn=min(XX);';
174         'ymx=max(YY);ymn=min(YY);';
175         'dx=xmx-xmn;if dx==0 then dx=max(xmx/2,1);end';
176         'xmn=xmn-dx/20;xmx=xmx+dx/20;';
177         'dy=ymx-ymn;if dy==0 then dy=max(ymx/2,1);end;';
178         'ymn=ymn-dy/20;ymx=ymx+dy/20;';
179         'xx2=orig(1)+sz(1)*((XX-xmn)/(xmx-xmn));';
180         'yy2=orig(2)+sz(2)*((YY-ymn)/(ymx-ymn));';
181         'xset(''color'',2)';
182         'xpoly(xx2,yy2,''lines'');']
183   x=standard_define([2 2],model,exprs,gr_i)
184 end
185 endfunction
186
187
188 function [rpar,ipar,ok] = poke_point(ixy,iparin,rparin)
189
190  [lhs,rhs]=argn(0)
191
192 //** get_click is already defined in "editi_curv"  
193 //in line definition of get_click
194 deff('[btn,xc,yc,win,Cmenu]=get_click(flag)',[
195 'if ~or(winsid() == curwin) then   Cmenu = ''Quit'';return,end,';
196 'if argn(2) == 1 then';
197 '  [btn, xc, yc, win, str] = xclick(flag);';
198 'else';
199 '  [btn, xc, yc, win, str] = xclick();';
200 'end;'; 
201 'if btn == -1000 then';
202 '  if win == curwin then';
203 '    Cmenu = ''Quit'';';
204 '  else';
205 '    Cmenu = ''Open/Set'';';
206 '  end,';
207 '  return,';
208 'end';
209 'if btn == -2 then';
210 '  xc = 0;yc = 0;';
211 '  try '    // added to handle unwanted menu actions in french version
212 '    execstr(''Cmenu='' + part(str, 9:length(str) - 1));';
213 '    execstr(''Cmenu='' + Cmenu);';
214 '  catch'
215 '    Cmenu=[]'    
216 '  end '    
217 '  return,';
218 'end';
219 'Cmenu=[]'])
220  
221 ok = %f
222 if rhs==0 then ixy=[];end;
223
224 if size(xy,'c')<2 then 
225   xinfo(" No [y] is provided");
226   return
227 end
228
229 [xy] = cleandata(ixy)
230
231 N = size(xy,'r');
232
233 if rhs<=1 then
234   NOrder = 1;
235   PeridicOption = 0;
236   ipar = [N;NOrder;PeridicOption]
237   rpar = []
238 else if rhs==2 then  
239     NOrder = iparin(2);
240     PeridicOption = iparin(3);
241     ipar = iparin;
242     rpar = []; 
243 else if rhs==3 then  
244     NOrder = iparin(2);
245     PeridicOption = iparin(3);
246     ipar = iparin;
247     rpar = rparin    
248 end
249 end //** ???
250 end //** ???
251
252 Amp=[]; wp=[]; phase=[]; offset=[]; np1=[];
253 Sin_exprs = list(string(Amp),string(wp), string(phase),string(offset),string(np1));
254 sAmp=[]; sTp=[]; sdelay=[];
255 Sawt1_exprs = list(string(sAmp),string(sTp),string(sdelay));
256 sAmp2=[]; sTp2 = [];
257 Sawt2_exprs = list(string(sAmp2),string(sTp2));
258
259 Amp3=[];Tp3=[];Pw3=[];Pd3=[];Bias3=[];
260 Pulse_exprs=list(string(Amp3), string(Tp3),string(Pw3),string(Pd3),string(Bias3))
261
262 mean4=[];var4=[];seed4=[];sample4=[];np4=[];
263 random_n_exprs=list(string(mean4),string(var4), string(seed4),string(sample4),string(np4))
264
265 min5=[];max5=[];seed5=[];sample5=[];np5=[];
266 random_u_exprs=list(string(min5), string(max5), string(seed5),string(sample5),string(np5))
267
268 // bornes initiales du graphique
269 xmx = max(xy(:,1));xmn=min(xy(:,1)),xmn=max(xmn,0);
270 ymx = max(xy(:,2));ymn=min(xy(:,2));
271 dx = xmx-xmn; dy = ymx-ymn
272
273 if dx==0 then dx=max(xmx/2,1),end;
274 xmx = xmx+dx/50;
275
276 if dy==0 then dy=max(ymx/2,1),end;
277 ymn = ymn-dy/50; ymx = ymx+dy/50;
278
279 rect = [xmn,ymn;xmx,ymx];
280
281 // initial draw
282 f = scf(curwin);
283
284
285 menu_r = [];
286 menu_s = [];
287 menu_o = ['zero order','linear','order 2','not_a_knot','periodic','monotone','fast','clamped']
288 menu_d = ['Clear','Data Bounds','Load from text file','Save to text file','Load from Excel','Periodic signal']
289 menu_t=['sine','sawtooth1','sawtooth2','pulse','random normal','random uniform']
290 menu_e=['Help','Exit without save','Save/Exit']
291 MENU=['Autoscale','Spline','Data','Standards','Exit'];
292 menus = list(MENU,menu_s,menu_o,menu_d,menu_t,menu_e);
293
294 scam='menus(1)(1)'
295 w='menus(3)(';r=')';Orderm=w(ones(menu_o))+string(1:size(menu_o,'*'))+r(ones(menu_o))
296 w='menus(4)(';r=')';Datam=w(ones(menu_d))+string(1:size(menu_d,'*'))+r(ones(menu_d))
297 w='menus(5)(';r=')';Standm=w(ones(menu_t))+string(1:size(menu_t,'*'))+r(ones(menu_t))
298 w='menus(6)(';r=')';Exitm=w(ones(menu_e))+string(1:size(menu_e,'*'))+r(ones(menu_e))
299
300 execstr('Autoscale_'+string(curwin)+'=scam')
301 execstr('Spline_'+string(curwin)+'=Orderm')
302 execstr('Data_'+string(curwin)+'=Datam')
303 execstr('Standards_'+string(curwin)+'=Standm')
304 execstr('Exit_'+string(curwin)+'=Exitm')
305
306 addmenu(curwin,MENU(1))
307 addmenu(curwin,MENU(2),menu_o)
308 addmenu(curwin,MENU(3),menu_d)
309 addmenu(curwin,MENU(4),menu_t)
310 addmenu(curwin,MENU(5),menu_e)
311 //===================================================================
312
313 drawlater();
314 a = gca();
315 a.data_bounds  = rect;
316 a.axes_visible = 'on';
317 a.clip_state   = 'on';
318 xtitle( '', 'time', 'Output' ) ; 
319 a.title.font_size=2;
320 a.title.font_style=4;
321 a.title.foreground=2;
322
323 a.grid=[2 2];
324 xpolys(xy(:,1),xy(:,2),[-1]);   //children(2)
325 xpolys(xy(:,1),xy(:,2),[5]);    //children(1)
326 splines = a.children(1).children
327 points  = a.children(2).children
328 //---------------------------------------
329 [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
330 drawnow();
331 // -- boucle principale
332
333 lines(0);
334 while %t then //=================================================
335   N = size(xy,'r');
336   [btn,xc,yc,win,Cmenu] = get_click(); //** see 
337   if ((win>0) & (win<>curwin)) then
338     Cmenu="Mouse click is Offside!";
339   end
340   if Cmenu==[] then Cmenu='edit',end
341   if (Cmenu=='Exit') |(Cmenu=='Quit' ) then, ipar=[];rpar=[];ok=%f;return; end
342   //-------------------------------------------------------------------
343   if ((Cmenu=='zero order') | (Cmenu=='linear') | (Cmenu=='order 2')| ...
344       (Cmenu=='not_a_knot')| (Cmenu=='periodic')| (Cmenu=='monotone')| ...
345       (Cmenu=='fast')| (Cmenu=='clamped')) then
346         
347     select  Cmenu
348      case 'zero order' then
349       NOrder=0;
350      case 'linear' then
351       NOrder=1;
352      case 'order 2' then
353       NOrder=2;
354      case 'not_a_knot' then
355       NOrder=3;
356      case 'periodic' then
357       NOrder=4;
358      case 'monotone' then
359       NOrder=5;
360      case 'fast' then
361       NOrder=6;
362      case 'clamped' then
363       NOrder=7;
364     end
365     ipar(2)=NOrder;
366    [rpar,ipar]=AutoScale(a,xy,ipar,rpar)  
367   end
368   //-------------------------------------------------------------------  
369   select Cmenu
370    case 'Data Bounds' then
371       rectx=findrect(a);
372       [mok, xmn1, xmx1, ymn1, ymx1] = scicos_getvalue('Enter new bounds',['xmin';'xmax'; 'ymin';'ymax'], ..
373           list('vec', 1,'vec', 1,'vec', 1,'vec', 1), string(rectx(:)))
374       //drawlater();
375       if mok then
376         if (xmn1 > xmx1 | ymn1 > ymx1) then
377           xinfo('Incorrect bounds')
378           mok=%f;
379         end
380         if xmn1<0 then
381           xinfo('X should be positive')
382           mok=%f;
383         end
384         if mok then
385           a.data_bounds=[xmn1, ymn1; xmx1, ymx1];
386         end
387       end
388       //drawnow();//show_pixmap(); 
389     //-------------------------------------------------------------------  
390    case 'Autoscale' then 
391     [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
392       //-------------------------------------------------------------------  
393    case 'Periodic signal' then 
394     if PeridicOption==1 then, ans0='y',else, ans0='n',end;
395     [mok,myans]=scicos_getvalue('Generating periodic signal',['y/n'],list('str',1),list(ans0));
396     if ((myans=='y')|(myans=='Y')) then,PeridicOption=1,else,PeridicOption=0;end;
397     ipar(3)=PeridicOption;
398     [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
399     //-------------------------------------------------------------------
400    case 'sine' then 
401     [mok,Amp,wp,phase,offset,np1,Sin_exprs2]=scicos_getvalue(' Sine parameters', ...
402                                 ['Amplitude';'Frequency(rad/sec)'; ...
403                     'Phase(rad)';'Bias';'number of points'],list('vec',1,'vec',1,'vec',1, ...
404                                            'vec',1,'vec',1),Sin_exprs)
405     if np1< 2 then np1=2;end
406     if mok & wp>0  then
407       NOrder=3;
408       ipar(2)=NOrder;
409       phase=atan(tan(phase));
410       xt=linspace(0,%pi*2/wp,np1)';
411       yt=Amp*sin(wp*xt+phase)+offset;
412       xy=[xt,yt];       
413       [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
414       Sin_exprs=Sin_exprs2
415     end
416     //-------------------------------------------------------------------
417    case 'sawtooth1' then 
418     [mok,sAmp,sTp,sdelay,Sawt1_exprs2]=scicos_getvalue('Sawtooth signal parameters', ...
419                           ['Amplitude';'Period';'delay'], ...
420                           list('vec',1,'vec',1,'vec',1),Sawt1_exprs)   
421     if mok & sTp>0 then
422       NOrder=1;
423       ipar(2)=NOrder;
424       if sdelay<sTp then 
425         xt=[0;sdelay;sTp];
426         yt=[0;0;sAmp];
427       else
428         xt=[0];
429         yt=[0];
430       end             
431       xy=[xt,yt];       
432       [rpar,ipar]=AutoScale(a,xy,ipar,rpar);
433       Sawt1_exprs=Sawt1_exprs2
434     end
435     //-------------------------------------------------------------------
436    case 'sawtooth2' then     
437     [mok,sAmp2,sTp2,Sawt2_exprs2]=scicos_getvalue('Sawtooth signal parameters', ...
438                           ['Amplitude';'Period'],list('vec',1,'vec',1),Sawt2_exprs)    
439     if mok & sTp2>0 then
440       NOrder=1;
441       ipar(2)=NOrder;
442       xt=[0;sTp2];
443       yt=[sAmp2;-sAmp2];
444       xy=[xt,yt];       
445       [rpar,ipar]=AutoScale(a,xy,ipar,rpar);
446       Sawt2_exprs=Sawt2_exprs2
447     end
448     //-------------------------------------------------------------------
449    case 'pulse' then
450     [mok,Amp3,Tp3,Pw3,Pd3,Bias3,Pulse_exprs2] = scicos_getvalue('Square wave pulse signal', ...
451                                                         ['Amplitude';'Period (sec)';'Pulse width(% o"+...
452                                                          "f period)';'Phase delay (sec)';'Bias'],list('vec',1, ...
453                                                          'vec',1,'vec',1,'vec',1,'vec',1),Pulse_exprs);        
454     if mok & Tp3>0  then
455       NOrder=0;
456       ipar(2)=NOrder;
457       if (Pd3>0) then xt=0;yt=Bias3;else xt=[];yt=[]; end
458       //otherwise there would be double points at 0
459       if Pd3<Tp3 then 
460         if Pw3>0 then 
461           xt=[xt;Pd3; Pw3*Tp3/100+Pd3;Tp3];
462           yt=[yt;Amp3+Bias3;Bias3;Bias3];       
463         else
464           xt=[0;Tp3];yt=[Bias3;Bias3];            
465         end      
466       else
467         xt=[0;Tp3];yt=[Bias3;Bias3];              
468       end
469       
470       xy=[xt,yt];       
471       [rpar,ipar]=AutoScale(a,xy,ipar,rpar);
472       Pulse_exprs=Pulse_exprs2;
473     end
474      //-------------------------------------------------------------------
475    case 'random normal' then
476     [mok,mean4,var4,seed4,sample4,np4,random_n_exprs2]=scicos_getvalue('Normal (Gaussian) random signal', ...
477                                     ['Mean';'Variance';'Initial seed';'Sample time';'Number of points'],list('vec',1, ...
478                                                   'vec',1,'vec',1,'vec', ...
479                                                   1,'vec',1),random_n_exprs)        
480     if mok & sample4>0 then
481       NOrder=0;
482       ipar(2)=NOrder;      
483       rand('normal');  rand('seed',seed4);
484       xt=0:sample4:sample4*(np4-1);xt=xt(:);
485       yt=mean4+sqrt(var4)*rand(np4,1);
486       xy=[xt,yt];       
487       [rpar,ipar]=AutoScale(a,xy,ipar,rpar);
488       random_n_exprs2=random_n_exprs;
489     end
490      //-------------------------------------------------------------------
491    case 'random uniform' then
492     [mok,min5,max5,seed5,sample5,np5,random_u_exprs2]=scicos_getvalue('Uniform random signal', ...
493                                     ['Minimum';'Maximum';'Initial seed';'Sample time';'Number of points'],list('vec',1, ...
494                                                   'vec',1,'vec',1,'vec', ...
495                                                   1,'vec',1),random_u_exprs)        
496     if mok & sample5>0 then
497       NOrder=0;
498       ipar(2)=NOrder;      
499        rand('uniform'); rand('seed',seed5);
500       xt=0:sample5:sample5*(np5-1);xt=xt(:);
501       yt=min5+(max5-min5)*rand(np5,1);
502       xy=[xt,yt];       
503       [rpar,ipar]=AutoScale(a,xy,ipar,rpar);
504       random_u_exprs2=random_u_exprs;
505
506     end   
507     //-------------------------------------------------------------------
508    case 'Save/Exit' then
509     NOrder=ipar(2);
510     PeridicOption=ipar(3);
511
512     METHOD=getmethod(NOrder);
513     if (METHOD=='periodic') then // periodic spline
514       xy(N,2)=xy(1,2);
515     end
516     
517     if (METHOD=='order 2' | METHOD=='not_a_knot'|METHOD=='periodic' | METHOD=='monotone'| METHOD=='fast' | METHOD=='clamped') then 
518       rpar=[xy(:,1);xy(:,2);rpar];
519     else
520       if (METHOD=='zero order'|METHOD=='linear')
521         rpar=[xy(:,1);xy(:,2);]
522       end
523     end
524     
525     ok=%t
526     delete(f);
527     return 
528     //-------------------------------------------------------------------
529    case 'Exit without save' then 
530     ipar=[];
531     rpar=[];
532     ok=%f
533     delete(f);
534     return
535     //-------------------------------------------------------------------
536    case 'Clear' then    
537     xy=[0,0];
538     NOrder=0;
539     ipar(2)=NOrder;
540     [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
541     //----------------------------------------------------------------
542    case 'Edit text data NOT IN USE' then
543     //  editvar xy;
544     [mok,xt,yt]=scicos_getvalue('Enter x and y data',['x';'y'],list('vec',-1,'vec',-1),list(strcat(sci2exp(xy(:,1))),strcat(sci2exp(xy(:,2)))));
545     if mok then,    
546       xy=[xt,yt];
547       [xy]=cleandata(xy), 
548       [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
549     end
550     //---------------------------------------------------------------  
551    case 'Help' then
552     t1='Mouse-left click: adding a new point'
553     t2='Mouse-right click: remove a point'
554     t3='Mouse-left double click: edit a point''s coordinates'
555     t4='Mouse-left button press/drag/release: move a  point'
556     t5='Change the window size: ''Data'' menu -> ''Databounds'''
557     messagebox([t1;t2;t3;t4;t5],"modal","info");
558     //---------------------------------------------------------------  
559    case 'Load from Excel' then
560     [tok,xytt]=ReadExcel()
561     if tok then
562       xy=xytt;
563       NOrder=1
564       ipar(2)=NOrder;
565       [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
566     end
567    //---------------------------------------------------------------       
568     case 'Load from text file' then
569     [tok,xytt]=ReadFromFile()
570     if tok then
571       xy=xytt;
572       NOrder=1
573       ipar(2)=NOrder;
574       [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
575     end
576    //---------------------------------------------------------------     
577    case 'Save to text file' then    
578     [sok]=SaveToFile(xy)
579     //---------------------------------------------------------------       
580    case 'Replot' then
581     if xy<>[] then 
582       drawlater();
583       points.data=xy;
584       [rpar,ipar]=drawSplin(a,xy,ipar,rpar);
585       drawnow()
586     end
587     //----------------------------------------------------------
588    case 'edit' then 
589     HIT=%f
590     if N<>0 then
591       xt=xy(:,1);yt=xy(:,2);
592       dist=((xt-ones(N,1)*xc))^2+((yt-ones(N,1)*yc))^2
593       [dca,k]=min(dist);
594       rectx=a.data_bounds;
595       ex=abs(rectx(2,1)-rectx(1,1))/80;
596       ey=abs(rectx(2,2)-rectx(1,2))/80;
597       if (abs(xc-xt(k))<ex & abs(yc-yt(k))<ey) then 
598         HIT=%t
599       end
600     end
601     
602     //_________________________
603   //  if ~((NOrder==-1|NOrder==-2|NOrder==-3|NOrder==-4)) then
604       if (~HIT)&(btn==0 | btn==3) then    // add point
605         if (xc>=0) then 
606           if (xc==0) then 
607             zz=find(x==0);
608             xy(zz,:)=[];
609           end 
610           xy=[xy;xc,yc];
611           [xtt,k2]=gsort(xy(:,1),'r','i');xy=xy(k2,:)
612           drawlater();
613           points.data=xy;
614           [rpar,ipar]=drawSplin(a,xy,ipar,rpar);  
615           drawnow();
616         end           
617       end
618       
619       if (HIT)&(btn==2 | btn==5) then  //   remove point
620         if (xy(k,1)>0) |( xy(k,1)==0 & (size(find(xy(:,1)==0),'*')>1)) then 
621           xy(k,:)=[];
622         end
623         drawlater();
624          points.data = xy;
625          [rpar,ipar] = drawSplin(a,xy,ipar,rpar);  
626         drawnow(); 
627       end   
628
629       if (HIT)&(btn==0) then             // move point
630         [xy,rpar,ipar] = movept(a,xy,ipar,rpar,k)   
631       end
632             
633       if (HIT)&(btn==10) then             // change data:: double click
634         [mok,xt,yt]=scicos_getvalue("Enter new x and y",['x';'y'],...
635                              list('vec',1,'vec',1),list(sci2exp(xy(k,1)),sci2exp(xy(k,2))));
636         if mok then 
637           xy(k,:) = [xt,yt];
638           [xy] = cleandata(xy)
639           drawlater();
640            points.data=xy;
641            [rpar,ipar]=AutoScale(a,xy,ipar,rpar) 
642           drawnow()
643         end
644       end
645
646   //  end
647     //_________________________________
648    
649   end
650   //----------------------------------------------------------
651 end
652 endfunction
653 //========================================================================
654 function [orpar,oipar] = drawSplin(a,xy,iipar,irpar)
655   N=size(xy,'r');// new size of xy
656   x=xy(:,1);  y=xy(:,2);
657   points=a.children(2).children
658   splines=a.children(1).children
659   order=iipar(2);
660   periodicoption=iipar(3);
661   orpar=irpar;
662    
663   METHOD=getmethod(order);
664   
665   if periodicoption==1 then PERIODIC='periodic, T='+string(x(N)-x(1));
666   else PERIODIC='aperiodic';end  
667   a.title.text=[string(N)+' points,  '+'Method: '+METHOD+',  '+PERIODIC];
668
669   if (N==0) then, return; end
670   if (N==1) then, order=0; end
671 //  NP=50;// number of intermediate points between two data points 
672   [X,Y,orpar]=Do_Spline(N,order,x,y);
673   if (periodicoption==1) then 
674     X=[X;X($)];
675     Y=[Y;Y(1)];
676   else
677     xmx=max(points.data(:,1));  xmn=min(points.data(:,1));
678     XMX=max(0,xmx); XMN=max(0,xmn);
679     xmx1=max(a.x_ticks.locations)
680     XMX=max(XMX,xmx1);    
681     X=[X;XMX];
682     Y=[Y;Y($)];
683   end
684   //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
685   splines.data=[X,Y];    
686   oipar=[N;iipar(2);periodicoption]
687 endfunction
688 //=============================================================
689 function [xyt,orpar,oipar]=movept(a,xy,iipar,irpar,k)
690 //on bouge un point existant
691   points=a.children(2).children
692   splines=a.children(1).children  
693   oipar=iipar
694   orpar=irpar
695   order=iipar(2);
696   x=xy(:,1);  y=xy(:,2);  
697   
698   if (x(k)==0) then 
699     zz=find(x==0);
700     x(zz)=[];y(zz)=[];
701     ZERO_POINT=%t
702   else
703     x(k)=[];
704     y(k)=[]; 
705     ZERO_POINT=%f
706   end 
707
708   btn=-1
709
710   while ~(btn==3 | btn==0| btn==10| btn==-5)
711     rep=xgetmouse([%t %t]); xc=rep(1);yc=rep(2);btn=rep(3);
712     if (ZERO_POINT) then 
713       xc=0;
714     else
715       if (xc<=0) then 
716         zz=find(x==0);
717         x(zz)=[];y(zz)=[];
718         ZERO_POINT=%t;
719         xc=0;
720       end
721     end
722   
723     xt=[x;xc];
724     yt=[y;yc];
725     [xt,k2]=gsort(xt,'r','i');yt=yt(k2)
726     xyt=[xt,yt];
727     
728     drawlater();
729       points.data=xyt;    
730       [orpar,oipar]=drawSplin(a,xyt,oipar,orpar);  
731     drawnow()
732   end
733
734 endfunction
735
736 //==========================================================
737 function rectx = findrect(a) 
738   splines=a.children(1).children  
739   points=a.children(2).children
740
741   if (points.data==[]) then 
742     rectx=a.data_bounds;
743     return;
744   end    
745
746
747   ymx1=max(splines.data(:,2));  ymn1=min(splines.data(:,2))
748   
749   xmx=max(points.data(:,1));xmn=min(points.data(:,1));
750   ymx=max(points.data(:,2));ymn=min(points.data(:,2));
751
752   
753   XMX=max(0,xmx);            XMN=max(0,xmn);
754   YMX=max(ymx,ymx1);  YMN=min(ymn,ymn1);
755   
756   dx=XMX-XMN;dy=YMX-YMN
757   if dx==0 then dx=max(XMX/2,1),end;
758   XMX=XMX+dx/50
759   if dy==0 then dy=max(YMX/2,1),end;
760   YMN=YMN-dy/50;YMX=YMX+dy/50;  
761   rectx=[XMN,YMN;XMX,YMX];
762 endfunction
763
764 //============================================================
765 function [tok,xyo]=ReadExcel()
766   TA=['A';'B';'C';'D';'E';'F';'G';'H';'I';'J';'K';'L';'M';'N';'O';'P'; ...
767       'Q';'R';'S';'T';'U';'V';'W';'X';'Y';'Z';'a';'b';'c';'d';'e';'f'; ...
768       'g';'h';'i';'j';'k';'l';'m';'n';'o';'p';'q';'r';'s';'t';'u';'v'; ...
769       'w';'x';'y';'z'];
770   TN=['0','1','2','3','4','5','6','7','8','9'];
771   xyo=[];tok=%f;
772   while %t
773     [zok,filen,sheetN,xa,ya]=scicos_getvalue('Excel data file ',['Filename';'Sheet #"+...
774                     " ';'X[start:Stop]';'Y[start:stop]'],list('str',1, ...
775                                                   'vec',1,'str',1, ...
776                                                   'str',1), ...
777                                       list(['Classeur1.xls'],['1'],['C5:C25'],['D5:D25']));   
778     if ~zok then break,end
779
780     try
781         [fd,SST,Sheetnames,Sheetpos] = xls_open(filen);
782     catch
783         xinfo('Scicos canot find the excel file:'+filen);
784         break;
785     end 
786     try  
787     N=size(Sheetnames,'*');
788     if ((sheetN<=N) &(sheetN>0)) then 
789       [Value,TextInd] = xls_read(fd,Sheetpos(sheetN))
790       mclose(fd)
791     end
792     xa=strsubst(xa,' ',''); px=strindex(xa,':'); 
793     ya=strsubst(ya,' ',''); py=strindex(ya,':');
794     x1=part(xa,1:px-1); x2=part(xa,px+1:length(xa));
795     y1=part(ya,1:py-1); y2=part(ya,py+1:length(ya));
796
797     x1p=min(strindex(x1,TN));
798     if x1p==[] then, xinfo('Bad address in X:'+x1); break, end
799     x11=part(x1,1:x1p-1);x12=part(x1,x1p:length(x1));
800     
801     x2p=min(strindex(x2,TN));
802     if x2p==[] then, xinfo('Bad address in X:'+x2); break, end
803     x21=part(x2,1:x2p-1);x22=part(x2,x2p:length(x2));
804     
805     y1p=min(strindex(y1,TN));
806     if y1p==[] then, xinfo('Bad address in Y:'+y1); break, end
807     y11=part(y1,1:y1p-1);y12=part(y1,y1p:length(y1));
808
809     y2p=min(strindex(y2,TN));
810     if y2p==[] then, xinfo('Bad address in Y:'+y2); break, end
811     y21=part(y2,1:y2p-1);y22=part(y2,y2p:length(y2));
812        
813     // x11 x12: x21 x22
814     
815     lx11=length(x11);lx21=length(x21);
816     ly11=length(y11);ly21=length(y21)
817     xstC=0;for i=1:lx11,xstC=xstC+modulo(find(TA==part(x11,lx11-i+1)),26)*26^(i-1);end
818     xenC=0;for i=1:lx21,xenC=xenC+modulo(find(TA==part(x21,lx21-i+1)),26)*26^(i-1);end
819     ystC=0;for i=1:ly11,ystC=ystC+modulo(find(TA==part(y11,ly11-i+1)),26)*26^(i-1);end
820     yenC=0;for i=1:ly11,yenC=yenC+modulo(find(TA==part(y21,ly21-i+1)),26)*26^(i-1);end
821         
822     xstR=evstr(x12);
823     xenR=evstr(x22);
824     ystR=evstr(y12);
825     yenR=evstr(y22);
826     
827     [mv,nv]=size(Value)
828
829     if ~(xstR<=mv & xstR>0 & xenR<=mv & xenR>0&ystR<=mv & ystR>0&yenR<=mv&yenR>0 ) then 
830       xinfo('error in Row data addresses'); break
831     end
832     if ~(xstC<=nv & xstC>0 & xenC<=nv & xenC>0&ystC<=nv & ystC>0&yenC<=nv&yenC>0 ) then 
833       xinfo('error in Column data addresses'); break
834     end
835      
836     xo=Value(min(xstR,xenR):max(xstR,xenR),min(xstC,xenC):max(xstC,xenC));
837     yo=Value(min(ystR,yenR):max(ystR,yenR),min(ystC,yenC):max(ystC,yenC));
838     [nx,mx]=size(xo);// adjusting the x and y size
839     [ny,my]=size(yo);
840     N=min(nx,ny);
841     xo=xo(1:N,:);
842     yo=yo(1:N,:);
843     
844     xyo=[xo,yo];
845     [xyo]=cleandata(xyo)
846     
847     tok=%t;break,
848     catch
849          xinfo(' Scicos cannot read your Excel file, please verify"+...
850                    " the parameters ');          
851          break
852     end  
853   end
854   
855 endfunction
856 //---------------------------------------------------------------
857 function [xyo]=cleandata(xye)
858   xe=xye(:,1)
859   ye=xye(:,2)
860   
861   [nx,mx]=size(xe);// adjusting the x and y size
862   [ny,my]=size(ye);
863   N=min(nx,ny);
864   xe=xe(1:N,:);
865   ye=ye(1:N,:);
866
867   // checking for NULL data
868   for i=1:N
869     if (xe(i)<>xe(i)) then 
870       xinfo('x contains no data:x('+string(i)+')'); 
871       return;
872     end
873     if (ye(i)<>ye(i)) then 
874       xinfo('Y contains no data:y('+string(i)+')'); 
875       return;
876     end      
877   end
878   zz=find(xe<0);xe(zz)=[];ye(zz)=[]
879   if (find(xe==0)==[]) then // add zero point
880     xe($+1)=0;
881     ye($+1)=0;
882   end
883   
884   [xo,k2]=gsort(xe,'r','i');
885   yo=ye(k2)    
886  
887   xyo=[xo,yo];
888 endfunction
889 //---------------------------------------------------------------
890 function  [orpar,oipar] = AutoScale(a,xy,inipar,inrpar)   
891   drawlater();    
892    oipar = inipar
893    orpar = inrpar
894    points = a.children(2).children
895    splines = a.children(1).children
896    points.data = xy;
897    splines.data = xy;
898    [orpar,oipar] = drawSplin(a,xy,oipar,orpar);
899    rectx=findrect(a);     
900    a.data_bounds = rectx;
901   drawnow()
902 endfunction
903 //============================
904 function METHOD = getmethod(order)
905   select order
906    case 0 then, METHOD='zero order'
907    case 1 then, METHOD='linear'
908    case 2 then, METHOD='order 2'
909    case 3 then, METHOD='not_a_knot'
910    case 4 then, METHOD='periodic'
911    case 5 then, METHOD='monotone'
912    case 6 then, METHOD='fast'
913    case 7 then, METHOD='clamped'
914   end
915 endfunction
916 //=======================================
917 function [sok,xye] = ReadFromFile()
918   xye=[];sok=%f;
919   while %t
920     [sok,filen,Cformat,Cx,Cy]=scicos_getvalue('Text data file ',['Filename';'Reading [C] f"+...
921                     "ormat';'Abscissa column';'Output"+...
922                     " column'],list('str',1,'str',1,'vec',1,'vec',1), ...
923                                       list(['mydatafile.dat'],['%g %g'],['1'],['2']));       
924     if ~sok then break,end
925     px=strindex(Cformat,'%');
926     NC=size(px,'*');    
927     if NC==[] then, xinfo("Bad format in reading data file");sok=%f;break;end
928     Lx=[];
929     try
930       fd=mopen(filen,'r');
931       Lx=mfscanf(-1,fd,Cformat);
932       mclose(fd);
933     catch
934       xinfo("Scicos canot open the data file: " + filen);
935       break;
936     end 
937
938     [nD,mD] = size(Lx);
939     if ((mD==0) | (nD==0)) then,  xinfo("No data read");sok=%f;break;end
940     if (mD<>NC) then, xinfo('Bad format');sok=%f;break;end
941     
942     xe=Lx(:,Cx);ye=Lx(:,Cy);
943     xye=[xe,ye];
944     [xye]=cleandata(xye)
945     sok=%t;break,
946  end 
947 endfunction
948 //=======================================
949 function [sok]=SaveToFile(xye)
950   xe=xye(:,1)
951   ye=xye(:,2)
952   sok=%f;
953   while %t
954     [sok,filen,Cformat]=scicos_getvalue('Text data file ',['Filename';'Writing [C] f"+...
955                     "ormat'],list('str',1,'str',1), ...
956                                       list(['mydatafile.dat'],['%g %g']));       
957     if ~sok then break,end
958     px=strindex(Cformat,'%');
959     NC=size(px,'*');    
960     if NC<>2 then, xinfo("Bad format in writing data file");sok=%f;break;end
961
962     Cformat=Cformat+'\n';
963     
964     try
965       fd=mopen(filen,'w');
966       mfprintf(fd,Cformat,xe,ye);
967       mclose(fd);
968     catch
969       xinfo('Scicos canot open the data file:'+filen);
970       break;
971     end 
972
973     sok=%t;break,
974  end 
975 endfunction
976 //=========================================================
977 function [X,Y,orpar]=Do_Spline(N,order,x,y)
978   X=[];Y=[];orpar=[];
979
980   METHOD=getmethod(order);
981
982   if (METHOD=='zero order') then 
983     X=x(1);Y=y(1);
984     for i=1:N-1
985       X=[X;x(i);x(i+1);x(i+1)];
986       Y=[Y;y(i);y(i);y(i+1)];
987     end
988     return
989   end    
990   //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
991   if (METHOD=='linear') then
992     X=[];
993     for i=1:N
994       X=[X;x(i)];
995       Y=[Y;y(i)];
996     end
997     return
998   end
999   //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1000   if (N<25) then NP=10;else
1001    if (N<50) then NP=5;else
1002     if (N<100) then NP=2;else
1003       if (N<200) then NP=1;else
1004         NP=0;end;end;end;
1005   end
1006   for i=1:N-1
1007     X=[X;linspace(x(i),x(i+1),NP+2)']; // pour tous sauf "linear" et "zero order"
1008   end
1009   //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1010   if (N>2) & (METHOD=='order 2') then
1011     Z=ORDER2(x,y);
1012     A=Z(1:N-1);
1013     B=Z(N:2*N-2);
1014     C=Z(2*N-1:3*N-3);
1015     
1016     for j=1:size(X,'*')
1017       for i=N-1:-1:1
1018         if X(j)>=x(i) then,break;end
1019       end
1020       Y(j)=A(i)*(X(j)-x(i))^2+B(i)*(X(j)-x(i))+C(i);
1021     end    
1022     orpar=matrix(Z,-1,1)   
1023   end  
1024    //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1025   if (METHOD=='not_a_knot') then
1026     try
1027     d = splin(x, y, METHOD);
1028     Y = interp(X, x, y, d);    
1029     orpar=d(:);
1030     catch
1031      xinfo('ERROR in SPLINE: '+METHOD)
1032     end
1033     
1034   end
1035    //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1036   if (METHOD=='periodic') then
1037     if y(1)<>y(N) then 
1038       y(N)=y(1)
1039     end
1040     try 
1041       d = splin(x, y,METHOD);
1042       Y = interp(X, x, y, d);  
1043       orpar=d(:);
1044     catch
1045     xinfo('ERROR in SPLINE: '+METHOD)
1046     end
1047   end
1048     //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1049   if (METHOD=='monotone' ) then
1050     try
1051       d = splin(x, y, METHOD);
1052       Y = interp(X, x, y, d);  
1053       orpar=d(:);
1054     catch
1055     xinfo('ERROR in SPLINE: '+METHOD)
1056     end
1057   
1058   end
1059    //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1060   if (METHOD=='fast') then
1061     try
1062       d = splin(x, y, METHOD);
1063       Y = interp(X, x, y, d);    
1064       orpar=d(:);
1065     catch
1066       xinfo('ERROR in SPLINE:  '+METHOD)    
1067     end  
1068   end  
1069    //!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
1070   if (METHOD=='clamped') then
1071     try
1072       d = splin(x, y, METHOD,[0;0]);
1073       Y = interp(X, x, y, d);    
1074       orpar=d(:);
1075     catch
1076       xinfo('ERROR in SPLINE: '+METHOD)    
1077     end
1078   end  
1079   
1080 endfunction
1081 //=================================================  
1082 function [Z]=ORDER2(x,y)
1083 N=size(x,'*')-1;
1084 A=zeros(3*N-1,N*3);
1085 B=zeros(3*N-1,1);
1086 for i=1:N
1087    j=3*(i-1)+1;
1088    A(j,i+2*N)=1;
1089    B(j)=y(i);
1090    A(j+1,i)=(x(i+1)-x(i))^2;
1091    A(j+1,i+N)=x(i+1)-x(i);
1092    A(j+1,i+2*N)=1;
1093    B(j+1)=y(i+1);
1094 end
1095
1096 for i=1:N-1
1097    j=3*(i-1)+1;
1098    A(j+2,i)=2*(x(i+1)-x(i));
1099    A(j+2,i+N)=1;   
1100    A(j+2,i+N+1)=-1;
1101 end
1102
1103 Q=zeros(3*N,3*N);
1104 for i=1:N
1105   Q(i,i)=4*(x(i+1)-x(i))^2
1106   Q(i,i+N)=2*(x(i+1)-x(i))
1107   Q(i+N,i)=2*(x(i+1)-x(i))
1108   Q(i+N,i+N)=1;
1109 end
1110
1111 At=[Q,A';A,zeros(3*N-1,3*N-1)]
1112 Bt=[zeros(3*N,1);B]
1113 Zt=At\Bt;
1114 Z=Zt(1:3*N,1)
1115 endfunction
1116 //===================================================