function [interv,values,probs]=lloyd_max(N)
%function y=lloyd_max(N)
%Iteratively obtains the optimum (Lloyd-Max)
%quantizer for a normally distributed source
%(unit variance) and outputs the optimum quantization levels
%in vector "interv" and the optimum
%quantization values in vector "values".
%N is the number of bits of quantization.
%It makes use of the fact that the density is symmetric
%around zero and only quantizes the positive axis. For the
%negative axis, the quantization intervals and values are
%the negatives of the ones obtained.
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Program written for EE 455 by C. Georghiades, Sept. 20, 1998%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%Compute the number of quantization levels, L:
L=2^N;
%V below is just an initial range and it can be set to any reasonable value.
V=2030;
%Define initial quantization intervals (uniformly distributed)
%teste
a=linspace(0,V,L/2);
a=[a,0];
%teste
%Start the iteration for the Lloyd-Max optimum quantizer:
%Use enough the norm of the error for xhat to determine convergence
xhat_old=zeros(1,L/2);
error=1;
while(error>1.0e-6)
%Compute the conditional mean estimates, xhat, corresponding to a
%(it calls function gauss_distr() which returns the distribution
%function of a normal variable):
b=gauss_distr(a);
c=b(2:length(b))-b(1:length(b)-1);
cm=(exp(-a(1:length(a)-1).^2/2)-exp(-a(2:length(a)).^2/2))/sqrt(2*pi);
xhat=cm./c;
error=norm(xhat-xhat_old);
xhat_old=xhat;
%Now update a using ai=(xhat_i+xhat_(i+1))/2:
a1=(xhat(2:length(xhat))+xhat(1:length(xhat)-1))/2;
%Append the boundaries of the quantization interval
a=[0,a1,Inf];
end
interv=a;
values=xhat;
probs=c;
评论0