clear all
%low pass filter design:
k  = 20000;
zeta = 1.01;
G = tf(k,[1 2*zeta*sqrt(k) k]);

W = logspace(-0.01,3,100)';
[M,P ]= bode(G,W);
M = M(:);
P = P(:);

figure
semilogx(W,M,'b.-');
%draw the cutoffs:
hold on;
plot(W,ones(size(W))*0.99,'b:',W,ones(size(W))*0.1,'b:')


%% bad design although it satisfies the stated specs:
k2  = 20000;
zeta2 = 0.5;
G2 = tf(k,[1 2*zeta2*sqrt(k2) k2]);

[M2,P2 ]= bode(G2,W);
M2 = M2(:);
P2 = P2(:);

figure
semilogx(W,M2,'b.-');
%draw the cutoffs:
hold on;
plot(W,ones(size(W))*0.99,'b:',W,ones(size(W))*0.1,'b:',W,ones(size(W))*1.01,'b:');


Ts = 1e-3;
time = [0:Ts:10]';

true_signal = 4 + sin(5*time +0.2);
sensor_data = true_signal + 0.5*randn(length(time),1);

y_f1 = lsim(G,sensor_data,time);
y_f2 = lsim(G2,sensor_data,time);

figure;
plot(time,true_signal,'b',time,y_f1,'r',time,y_f2,'g');

save('noisy_data.mat','sensor_data','time','true_signal','-v6','-mat');