lpcdemo.asv

自己写的关于lpc分析与合成的Demo程序

tim=20;
p=10;
[x,fs,bit]=auread('speechdemo.au');
N=fs*(tim/1000);
N=round(N);
avg=sum(x)/length(x);
x=x-avg;
w=window(@hamming,N);

disp(N);
    i=1;
    j=1;
   lt=round(length(x)/N)-1; 
   for i=1:lt
        y1=x(j:j+N-1);
        y=y1.*w;
        pr=0;
        max=abs(y(1));
        flag=0;
        count=0;
        sigseg=0;
        noiseseg=0;
       for s=1:N,
           pr=pr+abs(y(s));
           sigseg=sigseg+(y1(s).^2);
           if abs(y(s)) > max
               max=abs(y(s));
           end
           if y(s) < 0
               
               if flag == 0
                   count=count+1;
                   flag=1;
               end
           end
           if y(s) >=0
               
               if flag == 1
                   count=count+1;
                   flag=0;
               end
           end
                
    end
    m(i)=max;
    E(i)=pr;
    zcrs(i)=count;

     if E(i) < .85
      %  result(i)=0;%threshold for Ein silence is 0.85, and for voiced and unvoiced chk zero crossings=20
      %    end
   % if E(i) > .85
        if zcrs(i) > 20
            result(i)=1;%result=1 implies unvoiced
            f0(i)=0;
        else
            result(i)=2;%result=2 implies voiced
            f0(i)=pitchaut(N,fs,y1);
        end
        %----------
        
        %p=16;
        a(i,1:(p+1))=lpc(y1,p);
     X = fft(x,2^nextpow2(2*N-1));
Rx = ifft(abs(X.^2));
su=0;
for v=2:p
    su=su+a(v)*Rx((length(Rx)/2)+v);
end
gain=su+Rx((length(Rx)/2)+1);
    
 
%-------------
%synthesis

if(result(i)==2)
    nimp=floor((tim/1000)*f0(i));
    if nimp==0
        nimp=1;
    end
%impgen=linspace(1,N,nimp);
index=floor(N/nimp);
r=1;
impgen=zeros(1,N);

    impgen(1)=60;
if nimp > 1
  for u=2:nimp
    impgen(r+index)=60;
    r=r+index;
   end
end

else
    impgen=randn(1,N);
end
b=[1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0];
out=filter(b,a(i,1:(p+1)),impgen);
out=out.*gain;
output(j:j+N-1)=out;
j=j+N;
for s=1:length(out)
noiseseg=noiseseg+((y1(s)-real(out(s))).^2);
end
SNRSEG(i)=10*log10(sigseg/noiseseg);
end
switch tim
case 10
wavwrite(output,fs,bit,'lpc10_10.wav');
case 20
    wavwrite(output,fs,bit,'lpc10_20.wav');
case 30
    wavwrite(output,fs,bit,'lpc10_30.wav');
case 40
    wavwrite(output,fs,bit,'lpc10_40.wav');
otherwise
    disp('invalid input');
end
sig=0;
noise=0;
for g=1:length(output)
    sig=sig+(x(g).^2);
    noise=noise+((x(g)-real(output(g))).^2);
end
SNR=(sig/noise);
%AVSNRG=sum(real(SNRSEG));
su=0;
for g=1:length(SNRSEG)
    su=su+SNRSEG(g);
end
AVSNRG=su/length(SNRSEG);
disp('SNR=');disp(real(SNR));
disp('AVSNRG');disp(real(AVSNRG));
subplot(2,1,1);plot(x(1:length(output)));
subplot(2,1,2);plot(real(output));