www.gusucode.com > ​基于线性预测系数(LPC)的语音信号重构matlab源码程序 > LPC/lpc_analysis.m

    %%%%%%%%%%%%%%%%%%%%%%%%%%%%
%       LPC by whitejch    %
%       2006/02/02         %
%%%%%%%%%%%%%%%%%%%%%%%%%%%%
close all; clear all; clc;

news = fopen('hi.wav' ,'r');
x = fread(news , 'short')/2^15;
fclose(news);
%x=x(1:240);                                       % 30 msec (240 samples)
x=x(1:length(x)); 
LENGTH=length(x);                                 % length of x[n]
fs=8000;                                          % sampling rate 
n=0:1/fs:(LENGTH-1)/fs;   

% ---------------------------------------------------------------------- %

figure(1)

% subplot(221)
% original speech signal %
subplot(2,2,1),plot(n*1000, x),grid ,hold on
xlabel('Time[msec]');    ylabel('Amplitude');

% WIN_LEN=160;                                    % window length 160 (20msec)
N=256;                                            % N-point fft

% prediction error signal %
order=12;                                         % order
[a,g]=lpc(x,order);                               % predictor coefficients
est_x=filter([0 -a(2:end)],1,x);                  % Estimated signal
plot(n*1000,est_x,'r--'),hold off
title('Original Signal and LPC Estimated signal')
legend('Original Signal','LPC Estimate')

% ---------------------------------------------------------------------- %
% subplot(223)
error=x-est_x;                                    % Prediction error
subplot(2,2,3),plot(n*1000,error), grid;
xlabel('Time[msec]');    ylabel('Amplitude');
title('Prediction Error')

% ---------------------------------------------------------------------- %
% subplot(222)
% signal and LPC spectrum %
f=0:fs/N:fs/2-1;                                  % half frequency
X=abs(fft(x,N));                                  
X_LOG=20*log10(X);                                % LOG (db) of signal spectrum
% subplot(2,2,2),plot(f/1000,X_LOG(1:N/2)),hold on
% xlabel('frequency[kHz]'),ylabel('LOG(db)')

%dc=sum(x(1:240))/240;
EST_X=freqz(1,a,N/2);
EST_X_LOG=20*log10(abs(EST_X));                   % LPC spectrum
% plot(f/1000,EST_X_LOG(1:N/2),'r'),hold off
% title('Signal and LPC spectrum')

% ---------------------------------------------------------------------- %
% subplot(224)
% error spectrum %
ERROR=abs(fft(error,N));
ERROR_LOG=20*log10(ERROR);
subplot(2,2,4),plot(f/1000,ERROR_LOG(1:N/2)), grid
xlabel('frequency[kHz]'),ylabel('LOG(db)')
title('Error spectrum')

% ---------------------------------------------------------------------- %
% vocal tract function
% 2006/02/11

%figure(2)
subplot(222)

rr=xcorr(x);
%rr=rr(240:479);
rr=rr(length(x):end);
[a2,E]=levinson(rr,order);

H2=freqz(1,a2,N/2);                               
H2_dB=20*log10(abs(H2));                          
 
E_dB=10*log10(E);
plot(f/1000,X_LOG(1:N/2),f/1000,H2_dB+E_dB,'r'),grid;
xlabel('frequency[kHz]'),ylabel('dB');

wavwrite(est_x,'hello');