% function system()

%%
%% Adaptive equalizers excercise for the course
%% S-38.411 Signal Processing in Telecommunication I (2cr)
%%

%%%----------------------------------------------------------------------------
%%%				Initialization
%%%----------------------------------------------------------------------------
clear all;					% Clear all the variables.
K=500;						% Number of symbols.
I=1;						% Number of independent runs.
N=11;						% Number of filter coefficients.
MSE=0*ones(1,K);				% MSE for adaptive filter.
sigma2_n=0.001;					% Noise power.
hn=[0.2194 1.0000 0.2194];			% Channel impulse response.
my_max=1/(N*(sum(hn.^2)+sigma2_n));		% Maximum step size (Eq 5).
my=my_max/3;					% Step size parameter.
D=7;						% Delay factor.

%%%---------------------------------------------------------------------------
%%%				Transmitter
%%%---------------------------------------------------------------------------
x=sign(randn(1,K));				% Generate BPSK data.
ref_x=[0*ones(1,D),x];				% Generate reference signal.

%%%---------------------------------------------------------------------------
%%%				Channel
%%%---------------------------------------------------------------------------
for i=1:I,					% Run all the indepedent runs.
 ['Iteration:',mat2str(i)]			% Iteration index.
 channel_data=filter(hn,1,x);			% Distort the signal.
 n=sqrt(sigma2_n)*randn(1,K);			% Generate white noise sequence.
 channel_output=channel_data+n;			% AWGN channel.

 %%%---------------------------------------------------------------------------
 %%%				Equalizer
 %%%---------------------------------------------------------------------------
 [c,e,q,z]=lms(channel_output,ref_x,my,N);	% Linear LMS equalizer.
 MSE=MSE+abs(e).^2;				% Mean Square Error (MSE).
end

%%%---------------------------------------------------------------------------
%%%				Process the results
%%%---------------------------------------------------------------------------
if (I>1),

 % Average independent runs.
 MSE=MSE/I;

 % Plot MSE.
 title_string=['Parameters: N=',num2str(N),' D=',num2str(D),' my=',num2str(my,3)];
 semilogy(MSE,'g-');hold on;title(title_string);
 xlabel('TIME IN SAMPLES');ylabel('MSE');

 % Plot minimum MSE.
 semilogy(sigma2_n*ones(1,K),'b-.');

elseif (I==1)

 % Calculate the impulse responses.
 unit_imp=[1,0*ones(1,N-1)];
 ir_cha=filter(hn,1,unit_imp);
 ir_equ=filter(c,1,unit_imp);
 ir_sys=conv(ir_equ,ir_cha);
 max_y=max(max([ir_cha,ir_equ,ir_sys]))+0.1;
 min_y=min(min([ir_cha,ir_equ,ir_sys]))-0.1;range=(0:N-1);
 subplot(2,2,1);stem(range,ir_cha);axis([0 N-1 min_y max_y]);
 title('CHANNEL RESPONSE');xlabel('TIME IN SAMPLES');ylabel('AMPLITUDE');grid;
 subplot(2,2,2);stem(range,ir_equ);axis([0 N-1 min_y max_y]);
 title('EQUALIZER RESPONSE');xlabel('TIME IN SAMPLES');ylabel('AMPLITUDE');grid;
 subplot(2,2,3);stem((0:2*(N-1)),ir_sys);axis([0 2*(N-1) min_y max_y]);
 title('SYSTEM RESPONSE');xlabel('TIME IN SAMPLES');ylabel('AMPLITUDE');grid;
 pause;

 % Calculate the frequency responses.
 [freq_cha,freq_w]=freqz(ir_cha,1,512);
 [freq_equ,freq_w]=freqz(ir_equ,1,512);
 [freq_sys,freq_w]=freqz(ir_sys,1,512);freq_w=freq_w/(2*pi);
 max_y=max(max(abs([freq_cha,freq_equ,freq_sys])))+0.1;
 subplot(2,2,1);plot(freq_w,abs(freq_cha));
 title('FREQUENCY RESPONSE OF CHANNEL');xlabel('NORMALIZED FREQUENCY');
 ylabel('MAGNITUDE');axis([0 0.5 0 max_y]);grid;
 subplot(2,2,2);plot(freq_w,abs(freq_equ));
 title('FREQUENCY RESPONSE OF EQUALIZER');xlabel('NORMALIZED FREQUENCY');
 ylabel('MAGNITUDE');axis([0 0.5 0 max_y]);grid;
 subplot(2,2,3);plot(freq_w,abs(freq_sys));
 title('FREQUENCY RESPONSE OF CHANNEL+EQUALIZER');xlabel('NORMALIZED FREQUENCY');
 ylabel('MAGNITUDE');axis([0 0.5 0 max_y]);grid;
 pause;

 % Plot the signals: data, after channel, after equalizer, after slicer.
 max_y=max(max([x,channel_output,q,z]))+0.1;
 min_y=min(min([x,channel_output,q,z]))-0.1;range=(K-2*N:K);
 h=subplot(2,2,1);stem(range,x(range-D));
 title('BPSK SIGNAL');xlabel('TIME IN SAMPLES');
 ylabel('AMPLITUDE');axis([range(1) K min_y max_y]);grid;
 h=subplot(2,2,2);stem(range,channel_output(range-D));
 title('SIGNAL AFTER CHANNEL');xlabel('TIME IN SAMPLES');
 ylabel('AMPLITUDE');axis([range(1) K min_y max_y]);grid;
 h=subplot(2,2,3);stem(range,q(range));
 title('SIGNAL AFTER EQUALIZER');xlabel('TIME IN SAMPLES');
 ylabel('AMPLITUDE');axis([range(1) K min_y max_y]);grid;
 h=subplot(2,2,4);stem(range,z(range));
 title('SIGNAL AFTER SLICER');xlabel('TIME IN SAMPLES');
 ylabel('AMPLITUDE');axis([range(1) K min_y max_y]);grid;

end