Solution 8.8.3.mtr1x
The full MATLAB program used to nd the transfer function is:
%
%load the data
%
load mtr1step
t=mtr1step(1:500,1);;
y=mtr1step(1:500,2);;
plot(t,y)
pause
y=y(51:500);;
t=t(51:500);;
%y = y-y(1);;
y(1) = 0;;
t(1) = 0;;
plot(t,y)
pause
y(180:450) = mean(y(180:450));;
y(160:179) = mean(y(160:179));;
y(140:159) = mean(y(140:159));;
y(120:139) = mean(y(120:139));;
y(110:119) = mean(y(110:119));;
y(100:109) = mean(y(100:109));;
y(95:99) = mean(y(95:99));;
y(90:94) = mean(y(90:94));;
y(85:90) = mean(y(85:90));;
y(80:84) = mean(y(80:84));;
y(75:79) = mean(y(75:79));;
y(70:74) = mean(y(70:74));;
y(65:69) = mean(y(65:69));;
y(60:64) = mean(y(60:64));;
y(55:59) = mean(y(55:59));;
y(50:54) = mean(y(50:54));;
y(45:49) = mean(y(45:49));;
plot(t,y)
pause
%zero the summer
%
1
asum = 0.0;;
%
%Find the length of the array y(t) Note:this the step response. It has to be
%pulled into Matlab before this program is run.
%
top = size(y);;
top = top(1,1)
%
%Initialize counter.
%
k=1;;
%
%
%Find maximum value of step response.
%
K0 = max(y)
%
%
%Form function K0 - y(t). Note y(t) is the same as yu(t) in the write up.
%
ftemp = K0 - y;;
%
fsave1 = ftemp;;
%
%Reset counter
%
k=2;;
%
%Integrate K0 - y(t). and form the function y1(t).
%
while ( k <= top)
l= k-1;;
asum = asum + ( ( (ftemp(k)+ftemp(l))/2 ) *( t(k) - t(l) ) );;
t1(k) = t(k);;
f1(k) = asum;;
k=k+1;;
end
%
%Save data for plotting.
2
%
save time.dat t1 -ascii
save f1.dat f1 -ascii
save ftempa.dat ftemp -ascii
%
%SetK1 = to max(y1(t)), that is its final steady state value.
%
K1 = max(f1)
%
%
%form function: K1 - y1(t).
%
ftemp1 = K1 - f1;;
%
%Save function for later plotting
%
save ftempb.data ftemp1 -ascii
%
%reset summer
%
asum = 0.0;;
%
%Reset counter
%
k=2;;
%
%Integrate K1 - y1(t) to get y2(t) (called f2 in program)
%
while ( k <= top)
l= k-1;;
asum = asum + ( ( (ftemp1(k) + ftemp1(l) ) /2 ) *( t(k) - t(l) ) );;
f2(k) = asum;;
k=k+1;;
end
%
%SetK2 equal to maximum value of y2(t).
%
K2 = asum
%
%Save function y2(t) for later plotting.
3
%
save f2.data f2 -ascii
%
%Find poles of transfer function.
%
a1 = K1/K0
a2 = ( a1*K1 - K2 )/K0
cpol = [1, a1/a2 , 1/a2]
roots(cpol)
%
%Find gain of transfer function
%
Kplant = K0/a2
g=zpk([],[-3.6505 -8.2925],Kplant)
T=linspace(0,4.01,402);;
u=ones(402,1);;
[yhat,T] = lsim(g,u,T);;
%plot(t,y,'k-',T,yhat,'k--')
print -deps sr882c.eps
The rst step is to get the data into MATLAB. That is done by the MAT-
LAB statements:
load mtr1step
t=mtr1step(1:500,1);;
y=mtr1step(1:500,2);;
plot(t,y)
print -deps sr883mtr1xa.eps
The plot of the raw date is shown in Figure1. As can be seen, the plot
does not start at t =0and the initial voltage is not zero. The MATLAB
statements
y=y(51:500);;
t=t(51:500);;
%y = y-y(1);;
y(1) = 0;;
t(1) = 0;;
plot(t,y)
print -deps sr883mtr1xb.eps
produce the adjusted step respone shown in Figure 2 The rest of the MAT-
4
-0.1 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
0
0.5
1
1.5
2
2.5
3
3.5
Figure 1: Rawstep response data
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
0
0.5
1
1.5
2
2.5
3
3.5
Figure 2: Adjusted step response data
5
LAB program is designed to do the successiveintegrations necessary to nd
the twopoles. The MATLAB statements
y(180:450) = mean(y(180:450));;
y(160:179) = mean(y(160:179));;
y(140:159) = mean(y(140:159));;
y(120:139) = mean(y(120:139));;
y(110:119) = mean(y(110:119));;
y(100:109) = mean(y(100:109));;
y(95:99) = mean(y(95:99));;
y(90:94) = mean(y(90:94));;
y(85:90) = mean(y(85:90));;
y(80:84) = mean(y(80:84));;
y(75:79) = mean(y(75:79));;
y(70:74) = mean(y(70:74));;
y(65:69) = mean(y(65:69));;
y(60:64) = mean(y(60:64));;
y(55:59) = mean(y(55:59));;
y(50:54) = mean(y(50:54));;
y(45:49) = mean(y(45:49));;
plot(t,y)
print -deps sr883mtr1xc.eps
are a crude smoothing program designed to help the algorith converge. The
smoothed response is shown in Figure 3. Unfortunately the smoothing does
not work and the alogrithm does not converge to one pole in the righthalf
plane and one pole in the left half plane.
6
0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9
0
0.5
1
1.5
2
2.5
3
Figure 3: Smoth step response data
7