5a3e6604a75f91020592c520f7650c405e507f45
[scilab.git] / scilab / modules / cacsd / macros / bode.sci
1 // Scilab ( http://www.scilab.org/ ) - This file is part of Scilab
2 // Copyright (C) INRIA Serge Steer
3 // This file must be used under the terms of the CeCILL.
4 // This source file is licensed as described in the file COPYING, which
5 // you should have received as part of this distribution.  The terms
6 // are also available at    
7 // http://www.cecill.info/licences/Licence_CeCILL_V2-en.txt
8
9 function []=bode(varargin)
10   rhs=size(varargin)
11   if type(varargin($))==10 then
12     comments=varargin($),rhs=rhs-1;
13   else
14     comments=[];
15   end
16   fname="bode";//for error messages
17   fmax=[]
18   if or(typeof(varargin(1))==['state-space' 'rational']) then 
19     //sys,fmin,fmax [,pas] or sys,frq
20     refdim=1 //for error message
21     if rhs==1 then //sys
22       [frq,repf]=repfreq(varargin(1),1d-3,1d3)
23     elseif rhs==2 then //sys,frq
24       if size(varargin(2),2)<2 then
25         error(msprintf(_("%s: Wrong size for input argument #%d: A 1-by-n array expected with n>%d.\n"),..
26                        fname,1,1))
27       end
28       [frq,repf]=repfreq(varargin(1:rhs))
29     elseif or(rhs==(3:4)) then //sys,fmin,fmax [,pas]
30       [frq,repf]=repfreq(varargin(1:rhs))
31     else
32       error(msprintf(_("%s: Wrong number of input arguments: %d to %d expected.\n"),fname,1,5))
33     end
34     [phi,d]=phasemag(repf)
35     if rhs>=3 then fmax=varargin(3),end
36   elseif  type(varargin(1))==1 then 
37     //frq,db,phi [,comments] or frq, repf [,comments]
38     refdim=2
39     select rhs
40     case 2 then //frq,repf
41       frq=varargin(1);
42       if size(frq,2)<2 then
43         error(msprintf(_("%s: Wrong size for input argument #%d: A row vector with length>%d expected.\n"),..
44                        fname,1,1))
45       end
46       if size(frq,2)<>size(varargin(2),2) then
47         error(msprintf(_("%s: Incompatible input arguments #%d and #%d: Same column dimensions expected.\n"),..
48                          fname,1,2))
49       end
50       [phi,d]=phasemag(varargin(2))
51     case 3 then  //frq,db,phi
52       [frq,d,phi]=varargin(1:rhs)
53       if size(frq,2)<>size(d,2) then
54         error(msprintf(_("%s: Incompatible input arguments #%d and #%d: Same column dimensions expected.\n"),..
55                          fname,1,2))
56       end
57       if size(frq,2)<>size(phi,2) then
58         error(msprintf(_("%s: Incompatible input arguments #%d and #%d: Same column dimensions expected.\n"),..
59                          fname,1,3))
60       end
61     else 
62        error(msprintf(_("%s: Wrong number of input arguments: %d to %d expected.\n"),fname,2,4))
63     end
64   else
65     error(msprintf(_("%s: Wrong type for input argument #%d: Linear dynamical system or row vector of floats expected.\n"),fname,1))
66   end;
67    frq=frq';d=d',phi=phi'
68   [n,mn]=size(d)
69
70   if comments==[] then
71     hx=0.48
72   else
73     if size(comments,'*')<>mn then
74       error(msprintf(_("%s: Incompatible input arguments #%d and #%d: Same number of elements expected.\n"),...
75                      fname,refdim,rhs+1))
76     end
77     hx=0.43
78   end;
79
80   drawlater()
81   sciCurAxes=gca();
82   axes=sciCurAxes;
83   wrect=axes.axes_bounds;
84
85
86   //magnitude
87   axes.axes_bounds=[wrect(1)+0,wrect(2)+0,wrect(3)*1.0,wrect(4)*hx*0.95]
88   axes.data_bounds = [min(frq),min(d);max(frq),max(d)];
89   axes.log_flags = "lnn" ;
90   axes.grid=color('lightgrey')*ones(1,3);
91   axes.axes_visible="on";
92   if size(d,2)>1&size(frq,2)==1 then
93     xpolys(frq(:,ones(1,mn)),d,1:mn)
94   else
95     xpolys(frq,d,1:mn)
96   end
97   if fmax<>[]&max(frq)<fmax then
98       xpoly(max(frq)*[1;1],axes.y_ticks.locations([1 $]));e=gce();
99       e.foreground=5;
100   end
101   xtitle("",_("Frequency (Hz)"),_("Magnitude (Db)"));
102
103   //phase
104   axes=newaxes();
105   axes.axes_bounds=[wrect(1)+0,wrect(2)+wrect(4)*hx,wrect(3)*1.0,wrect(4)*hx*0.95];
106   axes.data_bounds = [mini(frq),mini(phi);maxi(frq),maxi(phi)];
107   axes.log_flags = "lnn" ;
108   axes.grid=color('lightgrey')*ones(1,3);
109   axes.axes_visible="on";
110   if size(phi,2)>1&size(frq,2)==1 then
111     xpolys(frq(:,ones(1,mn)),phi,1:mn)
112   else
113     xpolys(frq,phi,1:mn)
114   end
115   ephi=gce()
116
117   if fmax<>[]&max(frq)<fmax then
118       xpoly(max(frq)*[1;1],axes.y_ticks.locations([1 $]));e=gce();
119       e.foreground=5;
120   end
121   xtitle("",_("Frequency (Hz)"),_("Phase (°)"));
122   // create legend
123   if comments<>[] then
124     captions(ephi.children,comments,'lower_caption')
125   end
126   drawnow()
127   // return to the previous scale
128   set( "current_axes", sciCurAxes ) ;
129
130 endfunction