Solution 3.8.6.9
The rst step is to get the data into MATLAB. The followng commands do
that.
EDU>load tfid
EDU>size(tfid)
ans =
500 2
EDU>
the rst column is the time, and the second column is the step response.
Weusethe following commands to put the data in a form wecan plot the
step response versus time.
EDU>t = tfid(1:500,1);;
EDU>y = tfid(1:500,2);;
EDU>plot(t,y)
EDU>print -deps sr3869a.eps
EDU>
The time response is shown in Figure`1 As can be seen the time response,
whichwas taken o the oscilloscope shows negativetimebeforethe step is
applied and also the system does not start exactly at zero. The downward
trend just before t =0iswherethe switchwas thrown in the analog circuit
to initiate the step. So wehave some choices to make. We need to clean the
data up a bit to apply the identication technique described in this chapter.
First of all, we add 1.8 seconds to the t vector so that t(72) becomes 0.
Then wesety(72) equal to zero and t(72) equal to zero. Then wetruncate
the t and y vectors to the entries in locations 72 through vehundred. The
commands are:
t=t+1.8;;
t(72) = 0;;
y(72) = 0;;
t=t(72:500);;
y=y(72:500);;
plot(t,y);;
print -deps sr3869b.eps;;
1
-20 0 20 40 60 80 100
0
0.5
1
1.5
2
2.5
3
3.5
4
Figure 1: Rawstep response data
The adjusted response is shown in Figure 2. The rationale for the adjust-
ments is as follows. The response was generated bycreating a second order
comepensator on the Feedback and Control box(FAC). The system response
at steady state is almost exactly four in reponse to a 1 V input. Thus, all
wewanttodoisgetzero to be the pointatwhich the switchonthe FAC
was thrown to initiate the step and the response goes down close to zero.
This turns out to be approximately t = ;1:8s.Wearenowinaposition to
apply the identication techniques discussed in Chapter 3. As a rst step,
wecreate the function
y
1
= A;y(t)
and then
w(t)=ln(y(t));;
as shown bytheMATLAB dialogue:
EDU>max(y)
ans =
2
0 10 20 30 40 50 60 70 80 90
0
0.5
1
1.5
2
2.5
3
3.5
4
Figure 2: Adjusted step response
3.9800
EDU>y1 = 3.99-y;;
EDU>w = log(y1);;
EDU>plot(t,w)
EDU>print -deps lny1.eps
EDU>
Note that 3.98 was increased to 3.99 to avoid ln(0). The result is shown in
Figure 3. The slope is roughly
;
3:3
40
= ;0:0825
So there is a pole around s = ;0:1.
Our next goal is to nd B using the MATLAB dialogue:
EDU>B = (y-3.99)./exp(-0.082*t);;
EDU>plot(t(100:200),B(100:200))
EDU>mean(B(100:200))
3
0 10 20 30 40 50 60 70 80 90
-5
-4
-3
-2
-1
0
1
2
Figure 3: Plot of ln(1;y(t)
ans =
-4.4781
EDU>std(B(100:200))
ans =
0.4639
EDU>print -deps B.eps
EDU>
Note that wehave used the MATLAB command:
./
to do an elementbyelement division using two matrices. The plotofB
versus time is shown in Figure 4.
4
15 20 25 30 35 40
-6.5
-6
-5.5
-5
-4.5
-4
-3.5
-3
Figure 4: Comparison of adjusted response and y
1
(t)=1;e
;0:15t
Having found B we can nowattempt to nd p
2
.Todosoweform the
function
y
2
(t)=y(t);(A + Be
;p
1
t
= Ce
;p
2
t
:
In the presentcase
y
2
(t)=y(t);3:98+ 4:4781e
;0:0825t
:
We then take the natural logarithm of y
2
and plot it versus time. as shown
in Figure 5. Then the slope is
;
1:8;:5
2:5
= ;0:52:
Thus there is a pole at
s = ;0:52:
Then
C = ;(A + B)=;(3:98;4:4781) = 0:4981:
Then our estimate of y(t)is
^y(t)=3:98;4:4781e
;0:0825t
+0:4981e
;0:52t
:
5
0 2 4 6 8 10 12 14 16 18 20
-4.5
-4
-3.5
-3
-2.5
-2
-1.5
-1
-0.5
Figure 5: Plot of ln[y
2
(t)]
A comparison of y(t) and ^y(t)isshown in Figure 6. As can be seen the t
is very good.
The last taskistondK.Wehave
K
p
1
p
2
=3:98;;
whichinthe presentcasebecomes
K =3:980:08250:52 = 0:172:
Thus
G(s)=
0:172
(s +0:0:0825)(s+0:52)
:
6
0 10 20 30 40 50 60 70 80 90
0
0.5
1
1.5
2
2.5
3
3.5
4
Figure 6: Plot of ^y(t)versusy(t)
7