%Part 2
clear all
home
bwnm
bwms
plot(bwnm_dat(:,2),bwnm_dat(:,3));
[t,theta]=ginput(2)
wd1=(10*pi)/(t(2)-t(1));
wn1=wd1;
z1=((1/(t(2)-t(1)))*log(theta(1)/theta(2)))/wn1
ztest=1-z1^2

%Part 3
M=0.5;
d=0.04;
a=0.09;
Jm=4*((M*d^2)/8+(M*a^2))
plot(bwms_dat(:,2),bwms_dat(:,3));
[t2,theta2]=ginput(2)
wd2=(10*pi)/(t2(2)-t2(1));
wn2=wd2;
z2=((1/(t2(2)-t2(1)))*log(theta2(1)/theta2(2)))/wn2
J1=(Jm*wd2^2)/(wd1^2-wd2^2)
k=wn1^2*J1
D1=z1*wn1*2*J1
J2=J1+Jm
D2=z2*wn2*2*J2

%Part 4
t3=[0:0.001:4];
num1=[1 D1/J1];
den1=[1 D1/J1 k/J1];
num2=[1 D2/J2];
den2=[1 D2/J2 k/J2];
sys1=tf(num1,den1);
sys2=tf(num2,den2);
y1=impulse(sys1,t3);
y2=impulse(sys2,t3);

%Part 5
%plot(t3,y1*600)
%hold on
%plot(bwnm_dat(:,2),bwnm_dat(:,3))
%hold off
%pause
%hold on
plot(t3,y2*600)
hold on
plot(bwms_dat(:,2),bwms_dat(:,3))