%========= MMSE.M
% Trace les histogrammes du signal recu xk et 
% du signal apres egalisation ykwiener.
clear; N=5000; gc=[1 -1.4 0.8]; lg=length(gc);
ak=sign(randn(1,N)); % Suites des symboles
sk=filter(gc,1,ak);  % Signal emis
vsth=sqrt(gc*gc'); RSBdB=20;
%=====
sigma_b=vsth*10^(-RSBdB/20); bk=sigma_b*randn(1,N);
xk=sk+bk;            % Signal bruite recu
%===== MMSE (RIF a 50 coefficients)
LW=50; d=23;
rss=conv(gc,gc(lg:-1:1));
rsspos=[rss(lg:2*lg-1) zeros(1,LW-lg)];
Rxx=toeplitz(rsspos)+sigma_b*sigma_b*eye(LW);
ras=[zeros(1,d) gc(lg:-1:1) zeros(1,LW-lg-d)];
w=inv(Rxx)*ras';
ykwiener=filter(w,1,xk);
%====== RSB apres egalisation
roMMSE=max(conv(w,gc))/std(ykwiener);
%===== Affichage
points=50; subplot(211); hist(xk,points); grid
subplot(212); hist(ykwiener,points); grid