ex8_3b1.m

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

%%%%%%%%%%%%%%%% Example 8.3(b1) %%%%%%%%%%%%%%%%
%   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        %
%                September 2002                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   ----  Proportional-Integral-Derivative Control  ----
%
clear
disp('Example 8.3(b1)')

Ts = 0.1;			                % sampling period
Gp_s = tf(4,conv([2 1],[0.5 1]));   % continuous process
H_s = tf(1,[0.05 1]);		        % continuous sensor
Gp_z = c2d(Gp_s,Ts,'zoh');	        % discrete plant (no sensor)
HGp_s = H_s*Gp_s                    % continous plant
HGp_z = c2d(HGp_s,Ts,'zoh')         % discrete plant, using ZOH option

dtime  = 0:Ts:2;
KI = 0.5; 		
symbols = ['o' '+' '*' 'd' 'v'];    % different plot symbols
for KD = 0.4:0.05:0.6,
 figure
 results = [];
 Pgains = [];
 hold on
  for ii = 1:5,	                    % KP values are: 1.90,1.95,2.00,2.05,2.10
    KP = 1.85+0.05*ii;
    %-- write Gc(z) as sum of 3 LTI objects in TF form --
    Gc = KP*(tf(1,1,Ts) ...			% porportional
	+ tf([KI*Ts  0],[1 -1],Ts) ...	% integral
	+ tf(KD*[1 -1],[Ts 0],Ts)); 	% derivative

    % reference input - CL transfer fcn
    Tref_z = Gp_z*Gc/(1+HGp_z*Gc);
    Tref_z = minreal(Tref_z);       % remove pole-zero cancellations
    y = step(Tref_z,dtime);         % CL response to step input
    ref = dcgain(Tref_z);
    [Mo,tpeak,trise,tsettle,ess] = kstats(dtime,y,ref);
    ymax = y(find(dtime == tpeak));
    yss = y(length(dtime));
    results = [results; [KP KD trise tsettle  Mo]];
    Pgains = [Pgains KP];
    disp(['KP = ',num2str(KP)])
    plot(dtime,y,symbols(ii))
    pause(1)

  end  % for KP..
  hold off                          % turn hold off after last plot
  xlabel('Time (s)')
  ylabel('Step Response')
  title(['Example 8.3(b1) - Step Responses for Different ',...
        'Values for KP, where KI = ',num2str(KI),...
        ' and KD = ',num2str(KD)]) 
  grid on
  text('Units','norm','Pos',[0.2,0.15],'Str',['KI = ',num2str(KI)])
  text('Units','norm','Pos',[0.2,0.10],'Str',['KD = ',num2str(KD)])

  %------ show table of KP values and plot symbols -----
  text('Units','norm','Pos',[0.72,0.35],'Str','KP      symbol')
  text('Units','norm','Pos',[0.72,0.30],'Str','----------------------')
  text('Units','norm','Pos',[0.72,0.25],'Str','1.90      o')
  text('Units','norm','Pos',[0.72,0.20],'Str','1.95      +')
  text('Units','norm','Pos',[0.72,0.15],'Str','2.00      *')
  text('Units','norm','Pos',[0.72,0.10],'Str','2.05     diamond')
  text('Units','norm','Pos',[0.72,0.05],'Str','2.10     triangle')

  %----- print results table
  disp(['KI = ',num2str(KI),' and KD = ',num2str(KD)]) 
  disp(' ')
  disp('      KP       KD     t_rise      tsettle     M_o (pct)')
  disp('    ---------------------------------------------------')
  disp(results)
  disp('    ---------------------------------------------------')
  disp(' ')
  disp('******>'), pause
end  % for KD..
%%%%%%%%%%