% part 4 of the mini project: linear controller design:
% though no matlab script is asked for.
% Prabir Barooah, March 20, 2009

clear all
%conversion: 
mph2mpersec = 0.44704; %m/sec
mpersec2mph = 1/mph2mpersec;

m = 1300; %kg
b = 42926.2; %N/volt
a = 150.91; %30.17; % N/m^2 s^2 

v_0 = 60*mph2mpersec; %m/s
theta_0 = a*v_0^2/b; %volt

%more constants
p= 2*a*v_0/m;
k = b/m;

%PI controller design:
kp = 0.01 %something small
ki = (p+k*kp)^2/k/3 %just so that b^2-4ac<0

%make sure poles are complex
roots([ 1 (p+k*kp)   k*ki])
%at the same time we dont want the damping coefficient to be too small.


% construct the closed loop transfer function from r to y 
P = tf(k,[1 p]);
C = tf([kp ki],[1 0]);
H = feedback(series(P,C),1);
% if you try H = P*C/(1 + PC), it'll give you a transfer function that is a
% very very poor cousin of H. The reason is numerical (in) accuracy of the
% coefficients in the transfer functions P and C



%now, numerically test if the rise time and s.s. error indeed matches the specifications:
time = [0:0.001:10]';
y = step(H,time);
y_ss = 1; % steady state value of the step response of H

figure
plot(time,y)
%plot the 10% and 90% lines
hold on;
plot(time,y_ss*0.1*ones(length(time),1),':');
plot(time,y_ss*0.9*ones(length(time),1),':');
%to visually check the rise time

S = stepinfo(y,time,y_ss);