%	Figure 8.9
%	
%	BM
%	Dec 97


clear;
clf;
colordef(1,'black')

hopt = [ 6 4 ]';		%optimum weights
hint = [ 4 12 ]';		%initial guess at weights
hxr = [ 0 14 ]';			%range of tap values for contour plot
evr = [ 2 1 ];			%eigenvalues
theta = 30 ;			%orientation of eigenvectors (degrees)
sig_no = 10;			%signal/noise ratio (dB)
MM = 20;			%no. of iteration steps

%	set up optimum weights & correlation matrices
evr = evr .^2;
power = evr*ones(size(evr))';
evr = evr/power;
D = diag( evr );
theta = theta/180*pi;
V = zeros(2,2);
V(:,1) = [ cos(theta) -sin(theta) ]';
V(:,2) = [ sin(theta) cos(theta)]';
phiyy = V*D*V';					%autocorrelation matrix
sigpow = hopt'*phiyy*hopt;
phiyy = phiyy/sigpow;				%normalise o/p power
phiyx = phiyy*hopt;				%crosscorrelation vector
sig_no = 10^(sig_no/10);
nopow = 1/sig_no;				%noise power
sigx2 = 1+nopow;				%power in desired signal
mse_opt = 10*log10(sigx2 - hopt'*phiyx);	%optimum MSE in dB

%	mse surface
NN = 50;					
AA = zeros(NN,NN);
hmin = hxr(1);
hmax = hxr(2);
dh = (hmax-hmin)/(NN-1);
hx = hmin:dh:hmax;
hy = hx;
h = zeros(size(hopt));
for ii = 1:NN
	for jj = 1:NN
		h = [ hx(ii) hy(jj) ]';
		AA(ii,jj) = (sigx2 -2*h'*phiyx + h'*phiyy*h);
		end
	end

%	set up MSE surface plot
subplot(1,2,1)
contour(hx,hy,AA,20);
hold on
axis('square')
xlabel('tap 0');
ylabel('tap 1');

%	optimum weights
l0 = line('Xdata',hopt(2),'Ydata',hopt(1),'Marker','*','Erasemode','none','Color','w');

%	steepest descent
h = hint;
convg = zeros(2,MM);
convg(:,1) = h;
mse(1) = sigx2 -2*h'*phiyx + h'*phiyy*h;
mu = 1/3/trace(phiyy);
for ii = 2:MM
	grd = 2*(phiyy*h-phiyx);
	h = h - mu*grd;
	convg(:,ii) = h;
	mse(ii) = sigx2 -2*h'*phiyx + h'*phiyy*h;
end
mse = 10*log10(mse);

%	set up MSE convergence plot
subplot(1,2,2)
ii = 1;
m1 = line('Xdata',ii,'Ydata',mse(ii),'Marker','*','Color','y','Erasemode','xor');
m2 = line('Xdata',ii,'Ydata',mse(ii),'LineStyle','-','Color','r','Erasemode','none');
hold on
grid
axis ( [ 1 MM mse_opt max(mse) ]);
axis('square')
xlabel('iteration no');
ylabel('MSE (dB)');

subplot(1,2,1)
l1 = line('Xdata',convg(2,1),'Ydata',convg(1,1),'Marker','*','Color','y','Erasemode','xor');
l2 = line('Xdata',convg(2,1),'Ydata',convg(1,1),'LineStyle','-','Color','r','Erasemode','none');

%	ready to go
fprintf(1,'Figure 8.9: Press return\n')
pause

%	animation
for ii = 1:MM
	subplot(1,2,1)
	set(l1,'Xdata',convg(2,ii),'Ydata',convg(1,ii));
	if ii == 1,
		xx = [ convg(2,1) convg(2,1) ];
		yy = [ convg(1,1) convg(1,1) ];
	else
		xx = convg(2,ii-1:ii);
		yy = convg(1,ii-1:ii);
	end
	set(l2,'Xdata',xx,'Ydata',yy);
	
	subplot(1,2,2)
	set(m1,'Xdata',ii,'Ydata',mse(ii));
	if ii == 1,
		xx = [ ii ii ];
		yy = [ mse(ii) mse(ii) ];
	else
		xx = ii-1:ii;
		yy = mse(ii-1:ii);
	end
	set(m2,'Xdata',xx,'Ydata',yy);
	
	pause(1)

end