reinf6_5.m

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

%%%%%%%%%%% Reinforcement Problem 6.5 %%%%%%%%%%%
%   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                 %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%   ---- Chebyshev I bandpass filter ----
%
clear
disp('Reinforcement Problem 6.5')
Ts = 0.01;
%---- define system from Prob. 6.2
G = tf(0.0931*[1 0 -2 0 1],[1 0 1.0349 0 0.4293],Ts)

%-- first, draw Bode plot of filter to have as a guide ----
figure, bode(G)

%---- Part (a) --------------------
k1 = [0:1:100]';                      % generate time window
k1Ts = k1*Ts;
w1 = 4*pi                            % input frequency, rps
e1 = 2.0*sin(w1*k1Ts + 20*pi/180);    % input signal
y1 = lsim(G,e1);                     % simulate time response

[magw1,phw1] = bode(G,w1)                      % mag and phase at w1
y1ss = magw1*2*sin(w1*k1Ts + (20+phw1)*pi/180);  % eqn 6.3, but using sin
%--  magnitude and phase of steady-state y1 at w1
mag_y1_w1 = 2*magw1				% s-s magnitude 
ph_y1_w1 = 20+phw1				% s-s phase in degrees

figure
subplot(211)
plot(k1Ts,e1,'*');grid
xlabel('k1Ts (s)')
legend('e1')
title('Reinforcement Problem 6.5(a)')
subplot(212)
plot(k1Ts,y1,'*');grid
hold on
plot(k1Ts,y1ss,'o')
hold off
xlabel('k1Ts (s)')
legend('y1','y1ss',4)
%-- adjust axis so can see steady-state response
axis([0 k1Ts(end) -0.006 0.006]);  

%---- Part (b) --------------------
k2 = [0:1:30]';                        % generate appropriate time window
k2Ts = k2*Ts;
w2 = 60*pi                             % input frequency, rps
e2 = 2*cos(w2*k2Ts - 145*pi/180);      % input signal
y2 = lsim(G,e2);                       % simulate time response

[magw2,phw2] = bode(G,w2)              % mag and phase at w2 = 60 rps
y2ss = magw2*2*cos(w2*k2Ts + (-145+phw2)*pi/180);  % eqn 6.3
%--  magnitude and phase of steady-state y2 at w2
mag_y2_w2 = 2*magw2				% s-s magnitude 
ph_y2_w2 = -145+phw2				% s-s phase in degrees

figure
subplot(211)
plot(k2Ts,e2,'*');grid
hold on
plot(k2Ts,e2,':')
hold off
xlabel('k2Ts (s)')
legend('e2')
axis([0 0.30 -2.5 2.5])
title('Reinforcement Problem 6.5(b)')
subplot(212)
plot(k2Ts,y2,'*');grid
hold on
plot(k2Ts,y2,':')
plot(k2Ts,y2ss,'o')
hold off
axis([0 0.30 -2.5 2.5])
xlabel('k2Ts (s)')
legend('y2','y2ss')
%%%%%%%%%%