157 lines
3.8 KiB
Matlab
157 lines
3.8 KiB
Matlab
clc;
|
|
close all;
|
|
clear all;
|
|
|
|
% === Načítanie vstupnej stopy ===
|
|
[y, Fs] = audioread('sem2.wav');
|
|
t = (0:length(y)-1) / Fs;
|
|
|
|
figure; %1
|
|
plot(t, y);
|
|
xlabel('Time');
|
|
ylabel('Amplitude');
|
|
title('Waveform sem2.wav');
|
|
grid on;
|
|
|
|
[y1, Fs1] = audioread('flac.wav');
|
|
[y2, Fs2] = audioread('flac2.wav');
|
|
|
|
t1 = (0:length(y1)-1) / Fs1;
|
|
t2 = (0:length(y2)-1) / Fs2;
|
|
|
|
figure; %2
|
|
plot(t1, y1);
|
|
xlabel('Time [s]');
|
|
ylabel('Amplitude');
|
|
title('Waveform flac.wav');
|
|
grid on;
|
|
|
|
figure; %3
|
|
plot(t2, y2);
|
|
xlabel('Time [s]');
|
|
ylabel('Amplitude');
|
|
title('Waveform flac2.wav');
|
|
grid on;
|
|
|
|
% === FFT a spektrá ===
|
|
N1 = length(y1);
|
|
N2 = length(y2);
|
|
|
|
X1 = fft(y1);
|
|
X2 = fft(y2);
|
|
|
|
freq_shift1 = (-N1/2 : N1/2 - 1) * (Fs1 / N1);
|
|
freq_shift2 = (-N2/2 : N2/2 - 1) * (Fs2 / N2);
|
|
|
|
figure; %4
|
|
subplot(2,1,1);
|
|
plot(freq_shift1, real(fftshift(X1)), 'b');
|
|
title('Reálna časť posunutého spektra flac.wav');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Reálna hodnota');
|
|
grid on;
|
|
subplot(2,1,2);
|
|
plot(freq_shift1, imag(fftshift(X1)), 'r');
|
|
title('Imaginárna časť posunutého spektra flac.wav');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Imaginárna hodnota');
|
|
grid on;
|
|
|
|
figure; %5
|
|
subplot(2,1,1);
|
|
plot(freq_shift2, real(fftshift(X2)), 'b');
|
|
title('Reálna časť posunutého spektra flac2.wav');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Reálna hodnota');
|
|
grid on;
|
|
subplot(2,1,2);
|
|
plot(freq_shift2, imag(fftshift(X2)), 'r');
|
|
title('Imaginárna časť posunutého spektra flac2.wav');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Imaginárna hodnota');
|
|
grid on;
|
|
|
|
% === Nájdeme 5 dominujúcich frekvencií z flac2.wav ===
|
|
halfN2 = floor(N2/2);
|
|
realX2 = real(X2(1:halfN2));
|
|
imagX2 = imag(X2(1:halfN2));
|
|
freq2 = (0:halfN2-1)*(Fs2/N2);
|
|
|
|
[pks2r, locs2r] = findpeaks(abs(realX2), 'NPeaks', 3, 'SortStr', 'descend');
|
|
f0_2r = freq2(locs2r);
|
|
f0_2r(f0_2r < 2) = []; % odstráni DC
|
|
|
|
[pks2i, locs2i] = findpeaks(abs(imagX2), 'NPeaks', 2, 'SortStr', 'descend');
|
|
f0_2i = freq2(locs2i);
|
|
f0_2i(f0_2i < 2) = [];
|
|
|
|
disp('Najvýraznejšie frekvenčné zložky flac2.wav - reálna časť:');
|
|
disp(f0_2r');
|
|
disp('Najvýraznejšie frekvenčné zložky flac2.wav - imaginárna časť:');
|
|
disp(f0_2i');
|
|
|
|
% === IIR Notch Filter Design ==
|
|
f_all = unique([f0_2r(:); f0_2i(:)]);
|
|
if length(f_all) > 5
|
|
f_all = f_all(1:5);
|
|
end
|
|
|
|
r = 0.995; % stabilita
|
|
sos_all = [];
|
|
gain_all = 1;
|
|
|
|
for k = 1:length(f_all)
|
|
wo = f_all(k) / (Fs2 / 2);
|
|
bw = (1 - r);
|
|
[b_notch, a_notch] = iirnotch(wo, bw);
|
|
[sos, g] = tf2sos(b_notch, a_notch);
|
|
sos_all = [sos_all; sos]; %#ok<AGROW>
|
|
gain_all = gain_all * g;
|
|
end
|
|
|
|
% === Frekvenčná charakteristika filtra ===
|
|
figure;
|
|
freqz(sos_all, 1024, Fs2);
|
|
title('Frekvenčná charakteristika IIR Notch filtra (5 zložiek z flac2.wav)');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Zosilnenie [dB]');
|
|
grid on;
|
|
|
|
% === Filtrovanie flac2.wav ===
|
|
filtered_y2 = filtfilt(sos_all, gain_all, double(y2));
|
|
|
|
audiowrite('filtered_flac2.wav', filtered_y2, Fs2);
|
|
|
|
% === FFT po filtrovaní flac2 ===
|
|
X2_filt = fft(filtered_y2);
|
|
freq_shift2_filt = (-N2/2 : N2/2 - 1) * (Fs2 / N2);
|
|
|
|
figure;
|
|
subplot(2,1,1);
|
|
plot(freq_shift2_filt, real(fftshift(X2_filt)), 'b');
|
|
title('Reálna časť spektra po filtrovaní (filtered flac2.wav)');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Reálna hodnota');
|
|
grid on;
|
|
|
|
subplot(2,1,2);
|
|
plot(freq_shift2_filt, imag(fftshift(X2_filt)), 'r');
|
|
title('Imaginárna časť spektra po filtrovaní (filtered flac2.wav)');
|
|
xlabel('Frekvencia [Hz]');
|
|
ylabel('Imaginárna hodnota');
|
|
grid on;
|
|
|
|
% === Aplikácia filtra na pôvodný sem2.wav ===
|
|
filtered_sem2 = filtfilt(sos_all, gain_all, double(y));
|
|
audiowrite('filtered_sem2.wav', filtered_sem2, Fs);
|
|
|
|
% === Vizualizácia filtrovaného sem2 ===
|
|
figure;
|
|
plot(filtered_sem2);
|
|
title('Filtrovaný signál sem2.wav (filter z flac2.wav)');
|
|
xlabel('Vzorky');
|
|
ylabel('Amplitúda');
|
|
grid on;
|
|
|
|
disp('Filtrovanie hotové. Spočetné spektrá flac2 aj výsledný sem2 boli uložené.');
|