function[sqnr,a_quan,code]=u_pcm(a,n)
%u_pcm 一個序列的均勻PCM編碼
%[sqnr,a_quan,code]=u_pcm(a,n)
%a=輸入序列
%n=量化級數的數目（偶數）
%sqnr=輸入信號量化噪聲比（以dB為單位）
%a_quan=編碼前的量化輸出
%code=編碼后的輸出


amax=max(abs(a));
a_quan=a/amax;%歸一化
b_quan=a_quan;%賦初值
d=2/n;
q=d.*[0:n-1];
q=q-((n-1)/2)*d;%產生量化間隔
for i=1:n%下面進行量化
  a_quan(find((q(i)-d/2<=a_quan)&(a_quan<=q(i)+d/2)))=...
  q(i).*ones(1,length(find((q(i)-d/2<=a_quan)&(a_quan<=q(i)+d/2))));
  b_quan(find(a_quan==q(i))) =(i-1).*ones(1,length(find(a_quan==q(i))));
end
a_quan=a_quan*amax;%恢復原值
nu=ceil(log2(n));%以下時編碼,詳見課本
code=zeros(length(a),nu);
for i=1:length(a)
   for j=nu:-1:0
     if (fix(b_quan(i)/(2^j))==1)
       code(i,(nu-j))=1;
       b_quan(i)=b_quan(i)-2^j;
       end
   end
end
sqnr=20*log10(norm(a)/norm(a-a_quan));%求量化信噪比