% prj_1a.m Project #1 An extension of prj-1.m to (ii) allow input of parameters and 
%  (ii) repeated runs to compute a sample average of the mean-square error
%
% to observe the effect "simple moving average"
% on a noisy sequence, defined as
%	y(n) = (1/M)*sum[r(n-m), m=0,...,M-1], n=0,...,N-1
% where r(n)=x(r)+no(n)=signal + noise
%
% we encode the dsp statement into matlabese
%   by (i) setting n=nr-1 and (ii) m=mr-1  	
% the moving average now reads
% 	yr(nr) = (1/M)*sum[rr(nr-mr+1), mr=1,...,M], nr=1,...,N
% ***** NOTICE THE '+1' WHICH APPEARS IN THE ARGUMENT OF xr( ) !! ****************
% we can now implement this form in matlab
% we define x for more clarity
% this is encoded, then, into xr(nr)=1, nr in [1,16] and zero otherwise
%	
	clear,clc
	N=input('enter total sequence length N  ') 
	L=input('enter signal width L  ')
	D=input('enter signal delay D  ')			
	xr=[zeros(1,D),ones(1,L),zeros(1,N-L-D)];	% signal
	var=input('enter noise variance var  ')
	M=input('enter length of moving average  ')
%
	NT=input('enter number of sample runs NT  ')
for nt=1:NT
	no=randn(1,N);		% noise sequence
	rr=xr+sqrt(var)*no;	% noisy signal
% perform moving average
	for nr=1:N
		yr(nr)=0.0;
		for mr=1:M	% accumulate average
			if nr-mr+1>=1 & nr-mr+1<=N
				yr(nr)=yr(nr)+rr(nr-mr+1)/M;
			end
		end
	end
% compute sample average of MSEs
	MSE(nt)=sum((xr-yr).^2)/N;
end
SAVMSE=sum(MSE)/NT;
disp('sample averaged MSE is  ')
SAVMSE

