📄 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');
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -