added LPC/signal freq domain compare, cepstrum

This commit is contained in:
aj 2020-11-03 22:21:19 +00:00
parent 1828152664
commit 5eb36e7af8
10 changed files with 163 additions and 44 deletions

14
func/fft_.m Normal file
View File

@ -0,0 +1,14 @@
function [frequencies, values] = fft_(signal, sample_frequency)
L=length(signal);
Y = fft(signal);
P2 = abs(Y); % two-sided spectrum
% P2 = abs(Y/L); % two-sided spectrum
P1 = P2(1:L/2+1); % single-sided spectrum
P1(2:end-1) = 2*P1(2:end-1);
frequencies = sample_frequency*(0:(L/2))/L;
values = P1;
end

14
func/spectro.m Normal file
View File

@ -0,0 +1,14 @@
function spectro(signal, sample_frequency, windows, overlap_interval)
sample_overlap = (overlap_interval / 1000) * sample_frequency;
sample_size = size(signal);
%window_size = round(sample_size(1) / ((windows + 1)/2))
% Turn windows into window width in samples, take into account overlap
window_size = round((sample_size(1) + (windows + 1) * sample_overlap) / (windows+1));
spectrogram(signal, window_size, sample_overlap, [], sample_frequency, 'yaxis');
end

83
lpss.m Normal file
View File

@ -0,0 +1,83 @@
%% lpss.m
%%
%% Coursework script
close all;clear all;clc;
LPC_ORDER = 8;
DISPLAY_SAMPLES = 1000;
WINDOW_NUMBER = 10;
WINDOW_OVERLAP = 5; % ms
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% READ SIGNAL
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
[y, Fs] = audioread('samples/hood_m.wav');
L = length(y) % number of samples
DISPLAY_SAMPLES = min([DISPLAY_SAMPLES L]);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% LPC
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
a = aryule(y, LPC_ORDER) % signal, filter order
est_y = filter(0.02, a, y);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PREDICTION ERROR
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
e = y - est_y;
[acs, lags] = xcorr(e,'coeff');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% COMPARE TWO SIGNALS
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% x = 1:DISPLAY_SAMPLES;
% figure(1)
% plot(x, y(end-DISPLAY_SAMPLES+1:end), x, est_y(end-DISPLAY_SAMPLES+1:end), '--')
%plot(x, y(end-DISPLAY_SAMPLES+1:end))
%plot(x, est_y(end-DISPLAY_SAMPLES+1:end))
% grid
% xlabel('Sample Number')
% ylabel('Amplitude')
% legend('Original signal','LPC estimate')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PLOT AUTOCORRELATION
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% figure(2)
% plot(lags, acs)
% grid
% xlabel('Lags')
% ylabel('Normalized Autocorrelation')
%ylim([-0.1 1.1])
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PLOT FREQUENCY RESPONSE
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% CALCULATE FFT
[freq_dom_freqs, freq_dom_vals] = fft_(y, Fs);
% GET FILTER RESPONSE
[h, filter_freqs] = freqz(1, a, length(freq_dom_freqs), Fs);
figure(3)
plot(freq_dom_freqs, 20*log10(freq_dom_vals), 'r', filter_freqs, 20*log10(abs(h)), 'b')
% plot(w/pi, 20*log10(abs(h)))
grid
xlabel('Frequency (Hz)')
ylabel('Magnitude (dB)')
legend('Original Signal', 'LPC Filter')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PLOT ORIGINAL SPECTROGRAM
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% figure(4)
% spectro(y, Fs, WINDOW_NUMBER, WINDOW_OVERLAP);
% colormap bone
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% PLAY
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%sound(y, Fs);

View File

@ -3,12 +3,10 @@
%% Load wav and run through autocorr from econometrics toolbox
%% Used max range of m lags calculated from sample frequency
close all;
clear all;
close all;clear all;clc;
[y, Fs] = audioread('samples/hood_m.wav');
sample_size = size(y);
L = sample_size(1); % number of samples
L = length(y); % number of samples
f0 = 60; % low-pitched male speech
%f0 = 600; % children

14
lpss_cepstrum.m Normal file
View File

@ -0,0 +1,14 @@
%% lpss_cepstrum.m
%%
%% Load wav and calculate cepstrum
close all;clear all;clc;
% READ SIGNAL
[y, Fs] = audioread('samples/hood_m.wav');
t = (0:length(y) - 1);
c = rceps(y);
% c = cceps(y);
plot(t, c)
xlabel('Quefrency')

View File

@ -3,12 +3,10 @@
%% Load wav and plot in frequency domain using fft
%% Construct proper one-sided function for display
close all;
clear all;
close all;clear all;clc;
[y, Fs] = audioread('samples/hood_m.wav');
sample_size = size(y);
L = sample_size(1); % number of samples
L = length(y); % number of samples
% PLOT TIME DOMAIN
figure(1);

View File

@ -2,45 +2,48 @@
%%
%% Load wav and calculate LPC coeffs
close all;
clear all;
close all;clear all;clc;
ORDER = 40;
ORDER = 50;
DISPLAY_SAMPLES = 1000;
% READ SIGNAL
[y, Fs] = audioread('samples/hood_m.wav');
sample_size = size(y);
L = sample_size(1); % number of samples
L = length(y) % number of samples
DISPLAY_SAMPLES = min([DISPLAY_SAMPLES L]);
for ITER=1:ORDER
% CALCULATE FFT
[freq_dom_freqs, freq_dom_vals] = fft_(y, Fs);
for ITER=1:5:ORDER
% LPC
a = lpc(y,ITER); % signal, filter order
est_y = filter(0.02, a, y);
% PREDICTION ERROR
e = y - est_y;
[acs, lags] = xcorr(e,'coeff');
% COMPARE FREQ DOMAIN
[h, filter_freqs] = freqz(1, a, length(freq_dom_freqs), Fs);
% COMPARE TWO SIGNALS
x = 1:DISPLAY_SAMPLES;
figure(2)
plot(x, y(end-DISPLAY_SAMPLES+1:end), x, est_y(end-DISPLAY_SAMPLES+1:end), '--')
%plot(x, y(end-DISPLAY_SAMPLES+1:end))
%plot(x, est_y(end-DISPLAY_SAMPLES+1:end))
figure(1)
plot(freq_dom_freqs, 20*log10(freq_dom_vals), 'r', filter_freqs, 20*log10(abs(h)), 'b')
% plot(w/pi, 20*log10(abs(h)))
grid
xlabel('Sample Number')
ylabel('Amplitude')
legend('Original signal','LPC estimate')
xlabel('Frequency (Hz)')
ylabel('Magnitude (dB)')
legend('Original Signal', 'LPC Filter')
% PLOT AUTOCORRELATION
% COMPARE TWO SIGNALS TIME DOMAIN
% est_y = filter(0.02, a, y);
% x = 1:DISPLAY_SAMPLES;
% figure(2)
% plot(lags, acs)
% plot(x, y(end-DISPLAY_SAMPLES+1:end), x, est_y(end-DISPLAY_SAMPLES+1:end), '--')
% %plot(x, y(end-DISPLAY_SAMPLES+1:end))
% %plot(x, est_y(end-DISPLAY_SAMPLES+1:end))
% grid
% xlabel('Lags')
% ylabel('Normalized Autocorrelation')
% ylim([-0.1 1.1])
% xlabel('Sample Number')
% ylabel('Amplitude')
% legend('Original signal','LPC estimate')
pause(0.5)
end

View File

@ -2,8 +2,7 @@
%%
%% Load wav and play
close all;
clear all;
close all;clear all;clc;
[y, Fs] = audioread('samples/ee.wav');
size(y)

View File

@ -2,15 +2,13 @@
%%
%% Load wav and plot spectrogram
close all;
clear all;
close all;clear all;clc;
WINDOW_NUMBER = 20;
% matrix by audio channel and sample frequency
[y, Fs] = audioread('samples/hood_m.wav');
sample_size = size(y);
window_size = round(sample_size(1) / ((WINDOW_NUMBER + 1)/2));
window_size = round(length(y) / ((WINDOW_NUMBER + 1)/2));
freqRes = Fs / window_size
timeRes = window_size / Fs

View File

@ -2,18 +2,16 @@
%%
%% Load wav and plot spectrogram for iterating window sizes
close all;
clear all;
close all;clear all;clc;
%WINDOW_NUMBER = 20;
% matrix by audio channel and sample frequency
[y, Fs] = audioread('samples/hood_m.wav');
sample_size = size(y);
for WINDOW_NUMBER=1:50
window_size = round(sample_size(1) / ((WINDOW_NUMBER + 1)/2));
window_size = round(length(y) / ((WINDOW_NUMBER + 1)/2));
freqRes = Fs / window_size;
timeRes = window_size / Fs;