% demo of moving averages  
% J(Y)S

% resolution
N=1024;
t = (1:N);

% constant signal in noise
figure(1); set(1, 'Position', [  5 200 390 350]); title('noisy constant'); 
figure(2); set(2, 'Position', [400 200 390 350]); 
n=0.75;
x = 5*ones(size(t));
y = x + n*randn(size(t));
figure(1); plot ( t,y,'r', t,x,'y' ) ; axis([0 1024 0 10]);
'simple average'
z = sum(y)/N * ones(size(t));
figure(2); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);
pause;


% slowly varying signal in noise
figure(3); set(3, 'Position', [  5 200 390 350]); title('slowly varying signal in noise'); 
figure(4); set(4, 'Position', [400 200 390 350]); 
'simple moving average'
n=0.75;
x = 5 + 2*sin(t/256);
y = x + n*randn(size(t));
figure(3); plot ( t,y,'r', t,x,'y' ) ; axis([0 1024 0 10]);
L = 3;
f = ones(1,L)/L;
z = conv(y,f); z = z(floor(L/2):floor(L/2)+N-1);
figure(4); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);
L = 31;
f = ones(1,L)/L;
z = conv(y,f); z = z(floor(L/2):floor(L/2)+N-1);
figure(5); set(5, 'Position', [400 200 390 350]); 
figure(5); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);
L = 101;
f = ones(1,L)/L;
z = conv(y,f); z = z(floor(L/2):floor(L/2)+N-1);
figure(6); set(6, 'Position', [400 200 390 350]); 
figure(6); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);
pause;

% varying signal in noise
figure(7); set(7, 'Position', [  5 200 390 350]); title('signal in noise'); 
figure(8); set(8, 'Position', [400 200 390 350]); 
n=1.5;
x = 5 + 3*sin(t/4);
y = x + n*randn(size(t));
figure(7); plot ( t,y,'r', t,x,'y' ) ; axis([0 1024 0 10]);
'simple moving average can not work'
L = 101;
f = ones(1,L)/L;
z = conv(y,f); z = z(floor(L/2):floor(L/2)+N-1);
figure(8); set(8, 'Position', [5 200 790 350]); 
figure(8); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);
pause;
'1/4 1/2 1/4 filter is a good try'
L=3; f = [0.25 0.5 0.25];
z = conv(y,f); z = z(floor(L/2):floor(L/2)+N-1);
figure(9); set(9, 'Position', [ 5 200 790 350]); 
figure(9); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);
pause;
'properly designed long filter works well'
%L=101;f=fir1(100,0.1);
L=101; f=remez(100,[0,0.1,0.2,1],[1,1,0,0]);
z = conv(y,f); z = z(floor(L/2):floor(L/2)+N-1);
figure(10); set(10, 'Position', [ 5 200 790 350]);
figure(10); plot ( t,y,'r', t,z,'w' ) ; axis([0 1024 0 10]);