done_blk.m

离散控制系统设计的MATLAB 代码

%%%%%%%%%%%%%%%%% done_blk.m %%%%%%%%%%%%%%%%%%%%
%   Discrete-Time Control Problems using        %
%       MATLAB and the Control System Toolbox   %
%   by J.H. Chow, D.K. Frederick, & N.W. Chbat  %
%         Brooks/Cole Publishing Company        %
%                 November 2002                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

clear; close all
activity = -1;    % force activity menu to appear
while (activity < 1) | (activity > 3),
  disp(' '); disp...
  (' 1 <= Enter a discrete-time 1-block model (New, Unity-Feedback, etc.)')
  disp(' 2 <= Calculations (Poles, Zeros, Residues, etc.)')
  disp(' 3 <= Plots (Impulse, Step, Bode, Nyquist, etc.)'); disp(' ');
  activ = input('Select action, or enter 0 to quit: ','s');
  
if any(isletter(activ)) | str2num(activ) > 3,  
   disp(' '); disp('Need to enter a number between 0 and 3..');
   activity = -1;
else
    activity = str2num(activ);
end
    
if isempty(activity) | activity == 0, 
    disp('Completed single-block analysis'); 
    return 
end

if activity == 1,	% enter block model in either tf or zp form
model_form = -1;  % force model entry menu to appear
while (model_form < 1) | (model_form > 6),
  disp(' ');
  disp(' 1 <= Build new LTI object from TF data')
  disp(' 2 <= Build new LTI object from ZPK data')
  disp(' 3 <= Build new LTI object from SS data')
  disp(' 4 <= Display numerator and denominator values')
  disp(' 5 <= Display ZPK values')
  disp(' 6 <= Display SS values')
  disp(' 7 <= Create unity-feedback closed-loop model')
  disp(' 8 <= Activate keyboard'); disp(' ');
  model_form = input...
    ('Select model form or enter 0 for top-level menu: '); 

if (isempty(model_form) | model_form == 0),
  activity = -1;   % return to activity menu
  break   % terminate "while model_form..." loop
end    % of " if model_form == 0"

%+++++++++++++++++++ enter model in TF form ++++++++++++++++++++
  if model_form == 1,
    disp(' ')
  disp('   **** Interactive entry of transfer function using TF data ****')
  disp(' ')
% remove any stale zeros, poles, & gain from workspace
clear zer zer_str zer_last pol pol_str pol_last gain gain_str gain_last
clear Ts_last HOL
  disp...
  ('Enter coefficients of numerator and denominator polynomials of H(z)')
  disp('Default (initial) values are from Example 4.1')
  %----------------------------------------------------------------------
  num_last = [2 5.39 0.6012]*0.001; 
  den_last = [1 -2.01 1.17 -0.1272];
  Ts_last  = 0.1;
  num_str = vec2str(num_last);	% string representation of numerator
  den_str = vec2str(den_last);	% string representation of denominator
  Ts_str = num2str(Ts_last);
  cmd = 'n';	% 'n' ==> model is not OK
  while cmd == 'n',
    disp(' ')
    % allow user to enter new values for H(z)
    num = input(['Enter numerator - Default is "' num_str '" ==> ']);
    if isempty(num), num = num_last; end   %  CR ==> use old value
    den = input(['Enter denominator - Default is "' den_str '" ==> ']);
    if isempty(den), den = den_last; end   %  CR ==> use old value
    Ts = input(['Enter sampling period - Default is "' Ts_str '" ==> ']);
    if isempty(Ts), Ts = Ts_last; end % CR ==> use old value
    num_str = vec2str(num);  % 
    den_str = vec2str(den);
    Ts_str = num2str(Ts);
    disp('Numerator and denominator polynomials of TF form ')
    disp('    and sampling period are:')
    eval(['num = ' num_str]);
    eval(['den = ' den_str]);
    eval(['Ts = ' Ts_str]);
    if length(num) <= length(den), % proper systems 
       cmd = input('Is Model OK? (y|n) - Default:y - ','s');
       if isempty(cmd), 
          cmd = 'y'; disp('Model is OK, then proceed..'); 
       end
       num_last = num;  % save previous values in case want to reuse
       den_last = den;
       Ts_last = Ts;
       H = tf(num,den,Ts);
    else % check for improper systems
       disp(' '); disp('Improper System: deg(num) > deg(den). Re-enter.')
       disp(' '); cmd = 'n';
    end % if length(num)..
  end   % of "while" for "run another case"
  disp('Model has been entered as H, a discrete-time TF object')
  disp(' '); disp('H, in TF form, is: '); H, disp(' ');

  model_form = -1;   % force model_form menu to reappear
end  % of "if model_form ==1"

%++++++++++++++++ enter model in ZPK form +++++++++++++++++++++
if model_form ==2,
disp(' ')
disp('   **** Interactive entry of transfer function using ZPK data ****')
disp(' ')
disp('Enter H(z) in terms of zeros, poles, and gain ')
disp('Default values are from Example 4.1')
% remove any stale numerator & denominator from workspace
clear num num_str num_last pol pol_str pol_last Ts_last HOL
zer_last  = [-2.578 -0.1166]; 
pol_last  = [0.14 0.9350+0.1846i 0.9350-0.1846i];
gain_last = 0.002;
Ts_last   = 0.1;
zer_str = vec2str(zer_last);	% string representation of numerator
pol_str = vec2str(pol_last);	% string representation of denominator 
gain_str = num2str(gain_last);	% string representation of gain
Ts_str = num2str(Ts_last);	    %   "          "        "  Ts
cmd = 'n';	% 'n' ==> model not OK
while cmd == 'n',
  disp(' ')
  % allow user to enter new values for H(z)
  zer = input(['Enter zeros as column vector - Default: "' zer_str '" ==> ']);
  if isempty(zer), zer = zer_last; end   %  CR ==> use old value
  pol = input(['Enter poles as column vector - Default: "' pol_str '" ==> ']);
  if isempty(pol), pol = pol_last; end   %  CR ==> use old value
  gain = input(['Enter gain - Default: "' gain_str '" ==> ']);
  if isempty(gain), gain = gain_last; end   %  CR ==> use old value
  Ts = input(['Enter sampling period - Default: "' Ts_str '" ==> ']);
  if isempty(Ts), Ts = Ts_last; end		%  CR ==> use old value
  zer_str = vec2str(zer)  % 
  pol_str = vec2str(pol)
  gain_str = num2str(gain)
  Ts_str = num2str(Ts)
  disp('Zeros, poles, gain, and sampling period of ZPK form are:')
  eval(['zer = ' zer_str]);
  eval(['pol = ' pol_str]);
  eval(['gain = ' gain_str]);
  eval(['Ts = ' Ts_str]);
  cmd = input('Is Model OK? (y|n) - Default:y - ','s');
  if isempty(cmd), cmd = 'y'; disp('Model is OK, then proceed..'); end
  zer_last = zer;    % save previous values in case want to reuse
  pol_last = pol;
  gain_last = gain;
  Ts_last = Ts;
  H = zpk(zer,pol,gain,Ts);
end   % of "while cmd = n"
disp('Model has been entered as H, a discrete-time zero-pole-gain object')
disp(' '); disp('H, in ZPK form, is: '); H, disp(' ');
  model_form = -1;   % force model_form menu to reappear
end   % of "if model_form == 2"

%++++++++++++++++ enter model in SS form +++++++++++++++++++++
if model_form ==3,
disp(' ')
disp('   **** Interactive entry of elements of an LTI using SS data ****')
disp(' ')
disp('Enter H(z) in terms of matrices A, B, C, and D ')
disp('Default values are from Example 4.1')
% remove any stale numerator & denominator from workspace
clear num num_str num_last pol pol_str pol_last Ts_last HOL
A_last  = [0.98 0.095 0.003; -0.38 0.89 0.04; 0 0 0.14]; 
B_last  = [0.002; 0.05; 0.9];
C_last = [1 0 0];
D_last = 0;
Ts_last   = 0.1;
A_str  = mat2str(A_last);	% string representation of A-matrix
B_str  = vec2str(B_last);	% string representation of B-matrix 
C_str  = vec2str(C_last);	% string representation of C-matrix
D_str  = vec2str(D_last);	% string representation of D-matrix
Ts_str = num2str(Ts_last);  % string representation of  Ts
cmd = 'n';	% 'n' ==> model not OK
while cmd == 'n',
  disp(' ')
  % allow user to enter new values for H(z) - CR:carriage return
  A = input(['Enter A-matrix - Default: "' A_str '" ==> ']);
  if isempty(A), A = A_last; end    %  CR ==> use old value
  B = input(['Enter B-matrix - Default: "' B_str '" ==> ']);
  if isempty(B), B = B_last; end    %  CR ==> use old value
  C = input(['Enter C-matrix - Default: "' C_str '" ==> ']);
  if isempty(C), C = C_last; end    %  CR ==> use old value
  D = input(['Enter D-matrix - Default: "' D_str '" ==> ']);
  if isempty(D), D = D_last; end    %  CR ==> use old value
  Ts = input(['Enter sampling period - Default: "' Ts_str '" ==> ']);
  if isempty(Ts), Ts = Ts_last; end %  CR ==> use old value
  A_str = mat2str(A);  B_str = vec2str(B);
  C_str = vec2str(C);  D_str = vec2str(D);
  Ts_str = num2str(Ts);
  disp('A,B,C,D matrices and sampling period of H(z) in SS form are:')
  eval(['A = ' A_str]); eval(['B = ' B_str]);
  eval(['C = ' C_str]); eval(['D = ' D_str]);
  eval(['Ts = ' Ts_str]);
  cmd = input('Is Model OK? (y|n) - Default:y - ','s');
  if isempty(cmd), cmd = 'y'; disp('Model is OK, then proceed..'); end
  % save previous values in case want to reuse
  A_last = A; B_last = B; 
  C_last = C; D_last = D; 
  Ts_last = Ts;
  H = ss(A,B,C,D,Ts);
end   % of "while cmd = n"
disp('Model has been entered as H, a discrete-time state-space object')
disp(' '); disp('H, in SS form, is: '); H, disp(' ');
  model_form = -1;   % force model_form menu to reappear
end   % of "if model_form == 3"

%+++++++++++ Display numerator and denominator values ++++++++++++++
  if model_form == 4,
    if (exist('H')),
       disp('Numerator, denominator, and sampling period of H(z) are:')
       [num,den,Ts] = tfdata(H,'v')
       disp('********>'); pause
     else
       disp...
       ('Need model defined as LTI object to perform this opertaion')
       disp('*******>'); pause
    end
    model_form = -1;   % force model_form menu to reappear
  end    % of "model_form == 4"

%++++++++++ Display ZPK values ++++++++++++++
  if model_form == 5,
     if (exist('H')),
       disp('Zeros, poles, gain, and sampling period of H(z) are:')
       [zer,pol,gain,Ts] = zpkdata(H,'v')
       disp('********>'); pause
     else
       disp...
         ('Need model defined as LTI object to perform this operation')
       disp('*******>'); pause
     end    % of "if exist..."
     model_form = -1;   % force model_form menu to reappear
  end    % of "model_form == 5"
  
%++++++++++ Display SS values ++++++++++++++
  if model_form == 6,
     if (exist('H')),
       disp('A,B,C,D matrices and sampling period of H(z) are:')
       [A,B,C,D,Ts] = ssdata(H)
       disp('********>'); pause
     else
       disp...
         ('Need model defined as LTI object to perform this operation')
       disp('*******>'); pause
     end    % of "if exist..."
     model_form = -1;   % force model_form menu to reappear
  end    % of "model_form == 6"

%+++++++++++++++ unity feedback closed-loop system +++++++
 if model_form == 7,
   if (exist('H')),
       disp(' '); disp('Saving open-loop model as "HOL" ')
       disp('...and newly formed closed-loop system will now be "H"')
       disp(' '); 
       disp('Unity feedback closed-loop system, H(z), is:')
       HOL = H;
       H = feedback(HOL,1)
       disp('*******>');pause
     else
       disp(' ')
       disp('Need an open-loop model first to perform this operation..')
       disp('******>');pause
    end
    model_form = -1;   % force model_form menu to reappear
  end    % of "model_form == 7"
  
%+++++++++++++++++ activate keyboard ++++++++++++++++++++
if model_form == 8 ,
  disp('keyboard is active to allow:')
  disp('  * examining workspace')
  disp('  * entry of externally-defined model')
  disp('  * use of any Matlab command')
  disp('To resume, enter "return"')
  keyboard
  model_form = -1;   % return to model_form menu
end    % of "model_form == 8"
%-------
end   % of " while model_form....."
end   % of " if activity == 1"
%========== end of model creation ===========

%=========== calculations ===============
if activity == 2,
  if ~exist('H'),
      disp('Need to enter a model as LTI object first. Rerun.');
      break
  end % if ~exist('H')
  
  disp('Select calculations to be done for model')
  calc_type = -1;  % force calculation type menu to appear
  while (calc_type < 1) | (calc_type > 6),
    disp(' 1 <= Numerator, denominator, & sample period of H(z)')
    disp(' 2 <= Zeros, poles, gain, & sample period of H(z)')
    disp(' 3 <= Residues and poles of H(z)')
    disp(' 4 <= Display state-space matrix')
    disp(' 5 <= Damping ratio & natural freq of poles of H(z)')