clear all; close all; clc;

% spring damper parameters

M1 = 60/2.205; % converting Mass to kg
f1 = 1.1;
wn1 = 2*pi*f1;
K1 = M1*(wn1^2);
G1 = 0.01;
C1 = 2*((M1*K1)^0.5)*G1;

% pendulum damper parameters

m2_m1 = 0.05;
M2 = m2_m1*M1;
wn2 = 0.97*wn1;
K2 = M2*(wn2^2);
g = 9.81;
L = g/(wn2^2);
G2 = 0.15;
C2 = 2*((M2*K2)^0.5)*G2;

% TF for Q1

sys = tf([1 0 0],[M1 C1 K1]);

% Q2 (a)
figure(1);
impulse(sys);

% Q2 (b)
figure(2);
bode(sys);

% TF for Q3

s4 = M1*M2*L;
s3 = C2*L*(M1+M2) + C2*L*M2 - C1*M2*L;
s2 = M2*L*K1 + M2*g*(M1+M2) - C1*C2*L;
s1 = C2*L*K1 - M2*g*C1;
s = M2*g*K1;

sys2 = tf([M2*L C2*L M2*g 0 0],[s4 s3 s2 s1 s]);


% Q3 (b)
figure(3);
impulse(sys);
hold on;
impulse(sys2);
legend('Without PTMD','With PTMD');

% Q3 (b)
figure(4);
bode(sys);
hold on;
bode(sys2)
legend('Without PTMD','With PTMD');