% 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