clear all
clc
fs=5;
Ts=1/fs;
h=1/fs;   %時間步長
t=0:Ts:4095*Ts;
D=0.31;
a=1;
b=1;
s=0.3*sin(2*pi*0.01*t);
x1=s+sqrt(2*D)*randn(size(t));


y=fft(x1,4096); 
pyy=y.*conj(y)/4096;
ff=fs*(0:2048)/4096;


figure(1)
subplot(2,2,1)
plot(t,x1,'k')
title('Input signal sn')
xlabel('Time')
ylabel('Amplitude')
subplot(2,2,2)
plot(ff,pyy(1:2049),'k')
xlabel('Frequency ')
ylabel('Amplitude')
title('Spectrum of input signal sn')
axis([-0.1 0.1 0 1000])

x=sr(a,b,h,x1);

y=fft(x,4096); 
py=y.*conj(y)/4096;
ff=fs*(0:2048)/4096;


subplot(2,2,3)
plot(t,x,'k')
title('Output signal x')
xlabel('Time')
ylabel('Amplitude')
subplot(2,2,4)
plot(ff,py(1:2049),'k')
xlabel('Frequency')
ylabel('Amplitude')
title('Spectrum of output signal x')
axis([-0.1 0.1 0 1000])

% [tfr,t,f]=tfrsp(x');
% [tfr,rtfr,hat]=tfrrmsc(x');
% % colorbar;
% % figure(1)
% % contour(t,abs(f),tfr);
% figure(5)
% contour(t,abs(f)*10000,tfr);
% axis([0 1024 0 250])