sem2 done

semestralka2/analyza.m

@@ -0,0 +1,72 @@
+clc;
+close all;
+clear all;
+
+[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;
+
+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;
+set(gcf, 'Color', 'none'); set(gca, 'Color', 'none');
+
+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;
+set(gcf, 'Color', 'none'); set(gca, 'Color', 'none');

semestralka2/analyza_b.m

@@ -0,0 +1,114 @@
+clc;
+close all;
+clear all;
+
+[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;
+
+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;
+set(gcf, 'Color', 'none'); set(gca, 'Color', 'none');
+
+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;
+set(gcf, 'Color', 'none'); set(gca, 'Color', 'none');
+
+% === Najdenie špičiek len v prvej polovici spektra ===
+halfN1 = floor(N1/2);
+halfN2 = floor(N2/2);
+
+magX1 = abs(X1(1:halfN1));
+magX2 = abs(X2(1:halfN2));
+realX2 = real(X2(1:halfN2));
+imagX2 = imag(X2(1:halfN2));
+
+freq1 = (0:halfN1-1)*(Fs1/N1);
+freq2 = (0:halfN2-1)*(Fs2/N2);
+
+% ---- flac.wav (10 peaks from magnitude) ----
+[pks1, locs1] = findpeaks(magX1, 'NPeaks', 10, 'SortStr', 'descend');
+f0_1 = freq1(locs1);
+f0_1(f0_1 < 2) = []; % ignore DC and very low freq
+
+% ---- flac2.wav (4 peaks from real, 6 peaks from imag) ----
+[pks2r, locs2r] = findpeaks(abs(realX2), 'NPeaks', 4, 'SortStr', 'descend');
+f0_2r = freq2(locs2r);
+f0_2r(f0_2r < 2) = [];
+
+[pks2i, locs2i] = findpeaks(abs(imagX2), 'NPeaks', 6, 'SortStr', 'descend');
+f0_2i = freq2(locs2i);
+f0_2i(f0_2i < 2) = [];
+
+% ---- Výpis výsledkov ----
+disp('Najvýraznejšie frekvenčné zložky flac.wav (10, len 1/2 spektra):');
+disp(f0_1');
+
+disp('Najvýraznejšie frekvenčné zložky flac2.wav - reálna časť (4, len 1/2 spektra):');
+disp(f0_2r');
+
+disp('Najvýraznejšie frekvenčné zložky flac2.wav - imaginárna časť (6, len 1/2 spektra):');
+disp(f0_2i');
+
+
+
+
+% // NOTCH FILTER
+

semestralka2/analyza_c.m

@@ -0,0 +1,156 @@
+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é.');

semestralka2/analyza_d.m

@@ -0,0 +1,124 @@
+clc;
+close all;
+clear;
+
+%% === LOAD INPUTS ===
+[y, Fs] = audioread('sem2.wav');
+[y2, Fs2] = audioread('flac2.wav');
+
+if size(y,2) > 1, y = mean(y,2); end
+if size(y2,2) > 1, y2 = mean(y2,2); end
+
+t = (0:length(y)-1) / Fs;
+t2 = (0:length(y2)-1) / Fs2;
+
+%% === VISUALIZE ORIGINAL WAVEFORMS ===
+figure;
+subplot(2,1,1);
+plot(t, y);
+xlabel('Time [s]'); ylabel('Amplitude');
+title('sem2.wav – Original Signal');
+grid on;
+
+subplot(2,1,2);
+plot(t2, y2, 'r');
+xlabel('Time [s]'); ylabel('Amplitude');
+title('flac2.wav – Noise Sample');
+grid on;
+
+%% === FREQUENCY ANALYSIS (FFT of noise) ===
+N2 = length(y2);
+Y2 = fft(y2);
+f2 = (0:N2/2-1) * Fs2 / N2;
+
+magY2 = abs(Y2(1:N2/2));
+figure;
+plot(f2, 20*log10(magY2 + eps), 'r');
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+title('Noise Spectrum (flac2.wav)');
+xlim([0 Fs2/2]); grid on;
+
+%% === DETECT DOMINANT NOISE FREQUENCIES ===
+[pks, locs] = findpeaks(magY2, 'MinPeakHeight', max(magY2)*0.05, ...
+    'MinPeakDistance', round(10*N2/Fs2));
+f_peaks = f2(locs);
+f_peaks = f_peaks(f_peaks > 20 & f_peaks < Fs/2);
+f_peaks = sort(f_peaks(1:min(5,length(f_peaks))));
+fprintf('Detected noise frequencies [Hz]: %.1f ', f_peaks);
+fprintf('\n');
+
+%% === DESIGN CASCADE OF NOTCH FILTERS (z-plane zeros/poles) ===
+r = 0.995;
+sos_all = [];
+g_all = 1;
+
+for k = 1:length(f_peaks)
+    w0 = 2*pi*f_peaks(k)/Fs;
+    z = [exp(1j*w0) exp(-1j*w0)];
+    p = r*z;
+    b = real(poly(z));
+    a = real(poly(p));
+    b = b / (sum(b)/sum(a)); % normalize DC gain
+    [sos_i, g_i] = tf2sos(b, a);
+    sos_all = [sos_all; sos_i];
+    g_all = g_all * g_i;
+end
+
+%% === VISUALIZE FILTER FREQUENCY RESPONSE ===
+figure;
+freqz(sos_all, 2048, Fs);
+title('Cascade of Notch Filters (Pole-Zero Design)');
+
+%% === APPLY FILTER TO sem2.wav ===
+y_notch = filtfilt(sos_all, g_all, y);
+y_notch = y_notch / (max(abs(y_notch)) + eps) * 0.95;
+audiowrite('sem2_notch_filtered.wav', y_notch, Fs);
+fprintf('Saved: sem2_notch_filtered.wav\n');
+
+%% === ADDITIONAL SPEECH BANDPASS (100–8000 Hz) ===
+[z_bp, p_bp, k_bp] = butter(4, [100 8000]/(Fs/2), 'bandpass');
+[sos_bp, g_bp] = zp2sos(z_bp, p_bp, k_bp);
+y_final = filtfilt(sos_bp, g_bp, y_notch);
+y_final = y_final / (max(abs(y_final)) + eps) * 0.95;
+audiowrite('sem2_final.wav', y_final, Fs);
+fprintf('Saved: sem2_final.wav (Notch + Bandpass)\n');
+
+%% === TIME DOMAIN COMPARISON ===
+figure;
+subplot(3,1,1);
+plot(t, y, 'b'); title('Original sem2.wav');
+xlabel('Time [s]'); ylabel('Amplitude'); grid on;
+subplot(3,1,2);
+plot(t, y_notch, 'g'); title('After Notch Filtering');
+xlabel('Time [s]'); ylabel('Amplitude'); grid on;
+subplot(3,1,3);
+plot(t, y_final, 'r'); title('Final Output (Notch + Bandpass)');
+xlabel('Time [s]'); ylabel('Amplitude'); grid on;
+
+%% === FREQUENCY DOMAIN COMPARISON ===
+N = length(y);
+f = (0:N/2-1)*Fs/N;
+Y = fft(y);
+Y_notch = fft(y_notch);
+Y_final = fft(y_final);
+
+figure;
+subplot(3,1,1);
+plot(f, 20*log10(abs(Y(1:N/2))+eps), 'b');
+title('Original sem2.wav Spectrum');
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+xlim([0 2000]); grid on;
+subplot(3,1,2);
+plot(f, 20*log10(abs(Y_notch(1:N/2))+eps), 'g');
+title('After Notch Filtering Spectrum');
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+xlim([0 2000]); grid on;
+subplot(3,1,3);
+plot(f, 20*log10(abs(Y_final(1:N/2))+eps), 'r');
+title('After Notch + Bandpass Spectrum');
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+xlim([0 2000]); grid on;
+
+fprintf('\n Processing complete. Files saved:\n');
+fprintf('   sem2_notch_filtered.wav (Notch only)\n');
+fprintf('   sem2_final.wav (Notch + Bandpass)\n');

semestralka2/analyza_experimental.m

@@ -0,0 +1,265 @@
+clc;
+close all;
+clear all;
+
+%% === STAGE 0: Load Audio Files ===
+[y, Fs] = audioread('sem2.wav');
+[y1, Fs1] = audioread('flac.wav');
+[y2, Fs2] = audioread('flac2.wav');
+
+% Ensure mono
+if size(y,2) > 1, y = mean(y,2); end
+if size(y1,2) > 1, y1 = mean(y1,2); end
+if size(y2,2) > 1, y2 = mean(y2,2); end
+
+fprintf('Loaded audio files:\n');
+fprintf('  sem2.wav: %d samples, Fs = %d Hz\n', length(y), Fs);
+fprintf('  flac.wav: %d samples, Fs = %d Hz\n', length(y1), Fs1);
+fprintf('  flac2.wav: %d samples, Fs = %d Hz\n', length(y2), Fs2);
+
+%% === STAGE 1: IIR Notch Filter ===
+fprintf('\n=== STAGE 1: Notch Filtering ===\n');
+
+N2 = length(y2);
+X2 = fft(y2);
+halfN2 = floor(N2/2);
+freq2 = (0:halfN2-1) * (Fs2/N2);
+
+magX2 = abs(X2(1:halfN2));
+[pks, locs] = findpeaks(magX2, 'NPeaks', 10, 'SortStr', 'descend', 'MinPeakDistance', 5);
+f_peaks = freq2(locs);
+f_peaks = f_peaks(f_peaks > 10);
+if length(f_peaks) > 5
+    f_peaks = f_peaks(1:5);
+end
+
+fprintf('Detected tonal noise frequencies: ');
+fprintf('%.1f Hz, ', f_peaks);
+fprintf('\n');
+
+r = 0.995;
+sos_all = [];
+gain_all = 1;
+
+for k = 1:length(f_peaks)
+    wo = f_peaks(k) / (Fs2/2);
+    if wo > 0 && wo < 1
+        bw = 0.005;
+        [b_notch, a_notch] = iirnotch(wo, bw);
+        [sos, g] = tf2sos(b_notch, a_notch);
+        sos_all = [sos_all; sos];
+        gain_all = gain_all * g;
+    end
+end
+
+if ~isempty(sos_all)
+    y_notch = filtfilt(sos_all, gain_all, double(y));
+    y2_notch = filtfilt(sos_all, gain_all, double(y2));
+else
+    y_notch = double(y);
+    y2_notch = double(y2);
+end
+
+fprintf('Notch filtering complete.\n');
+
+%% === STAGE 2: Spectral Subtraction ===
+fprintf('\n=== STAGE 2: Spectral Subtraction ===\n');
+params.frameLen = 2048;
+params.hopLen = 512;
+params.nfft = 2048;
+params.alpha = 4.0;
+params.beta = 0.002;
+params.gamma = 1.0;
+[y_specsub, ~] = spectral_subtraction(y_notch, y2_notch, Fs, params);
+fprintf('Spectral subtraction complete.\n');
+
+%% === STAGE 3: MMSE-LSA Enhancement ===
+fprintf('\n=== STAGE 3: MMSE-LSA Enhancement ===\n');
+params_mmse.frameLen = 2048;
+params_mmse.hopLen = 512;
+params_mmse.nfft = 2048;
+params_mmse.alphaDD = 0.98;
+params_mmse.minGain = 0.1;
+params_mmse.noisePowSmooth = 0.7;
+[y_mmse, ~] = mmse_lsa_denoise(y_specsub, y2_notch, Fs, params_mmse);
+fprintf('MMSE-LSA enhancement complete.\n');
+
+%% === STAGE 4: Post-processing ===
+fprintf('\n=== STAGE 4: Post-processing ===\n');
+[b_hp, a_hp] = butter(2, 80/(Fs/2), 'high');
+y_final = filtfilt(b_hp, a_hp, y_mmse);
+y_final = y_final / (max(abs(y_final)) + eps) * 0.95;
+fprintf('Post-processing complete.\n');
+
+%% === Save outputs ===
+audiowrite('sem2_stage1_notch.wav', y_notch/max(abs(y_notch)+eps)*0.95, Fs);
+audiowrite('sem2_stage2_specsub.wav', y_specsub/max(abs(y_specsub)+eps)*0.95, Fs);
+audiowrite('sem2_stage3_mmse.wav', y_mmse/max(abs(y_mmse)+eps)*0.95, Fs);
+audiowrite('sem2_final_cleaned.wav', y_final, Fs);
+
+y2_cleaned = spectral_subtraction(y2_notch, y2_notch, Fs2, params);
+y2_cleaned = mmse_lsa_denoise(y2_cleaned, y2_notch, Fs2, params_mmse);
+y2_cleaned = y2_cleaned / (max(abs(y2_cleaned)) + eps) * 0.95;
+audiowrite('flac2_final_cleaned.wav', y2_cleaned, Fs2);
+fprintf('\n✅ All files saved!\n');
+
+%% === SIMPLE VISUALIZATION (same style as older system) ===
+t = (0:length(y)-1)/Fs;
+
+figure;
+subplot(3,1,1);
+plot(t, y); grid on;
+title('Original sem2.wav');
+xlabel('Time [s]'); ylabel('Amplitude');
+
+subplot(3,1,2);
+plot(t, y_notch); grid on;
+title('After Notch Filtering');
+xlabel('Time [s]'); ylabel('Amplitude');
+
+subplot(3,1,3);
+plot(t, y_final); grid on;
+title('Final Cleaned sem2\_final\_cleaned.wav');
+xlabel('Time [s]'); ylabel('Amplitude');
+
+figure;
+subplot(3,1,1);
+plot((0:length(y2)-1)/Fs2, y2); grid on;
+title('Original flac2.wav (Noise Sample)');
+xlabel('Time [s]'); ylabel('Amplitude');
+
+subplot(3,1,2);
+plot((0:length(y2_notch)-1)/Fs2, y2_notch); grid on;
+title('After Notch Filtering (Noise)');
+xlabel('Time [s]'); ylabel('Amplitude');
+
+subplot(3,1,3);
+plot((0:length(y2_cleaned)-1)/Fs2, y2_cleaned); grid on;
+title('Cleaned flac2.wav (Noise Reference)');
+xlabel('Time [s]'); ylabel('Amplitude');
+
+fprintf('\nProcessing summary:\n');
+fprintf(' - Stage 1 removed %d tonal frequencies\n', length(f_peaks));
+fprintf(' - Stage 2 applied spectral subtraction (α = %.1f)\n', params.alpha);
+fprintf(' - Stage 3 applied MMSE-LSA (α_DD = %.2f)\n', params_mmse.alphaDD);
+fprintf(' - Final result saved as sem2_final_cleaned.wav (best quality)\n');
+
+
+function [y_out, noise_psd] = spectral_subtraction(x, noise_ref, Fs, p)
+    x = x(:);
+    noise_ref = noise_ref(:);
+    frameLen = p.frameLen; hopLen = p.hopLen; nfft = p.nfft;
+    alpha = p.alpha; beta = p.beta; gamma = p.gamma;
+    win = hamming(frameLen, 'periodic');
+
+    [N_stft, ~, ~] = stft(noise_ref, Fs, 'Window', win, ...
+        'OverlapLength', frameLen-hopLen, 'FFTLength', nfft);
+    noise_psd = mean(abs(N_stft).^2, 2);
+    noise_psd = movmean(noise_psd, 5);
+
+    [X_stft, ~, ~] = stft(x, Fs, 'Window', win, ...
+        'OverlapLength', frameLen-hopLen, 'FFTLength', nfft);
+
+    [nFreq, nFrames] = size(X_stft);
+    Y_stft = zeros(size(X_stft));
+
+    for k = 1:nFrames
+        X_frame = X_stft(:, k);
+        X_mag = abs(X_frame);
+        X_phase = angle(X_frame);
+        X_pow = X_mag.^2;
+        snr_post = X_pow ./ (noise_psd + eps);
+        alpha_freq = alpha * (1 + 0.5 * (1 - min(snr_post, 4)/4));
+        Y_pow = X_pow - alpha_freq .* noise_psd;
+        Y_pow = max(Y_pow, beta * X_pow);
+        Y_mag = sqrt(max(Y_pow, 0));
+        Y_stft(:, k) = Y_mag .* exp(1j * X_phase);
+    end
+
+    y_out = istft(Y_stft, Fs, 'Window', win, ...
+                  'OverlapLength', frameLen-hopLen, ...
+                  'FFTLength', nfft, 'ConjugateSymmetric', true);
+    y_out = y_out(1:min(length(y_out), length(x)));
+    if length(y_out) < length(x)
+        y_out = [y_out; zeros(length(x)-length(y_out),1)];
+    end
+end
+
+
+function [y_out, G_all] = mmse_lsa_denoise(x, noise_ref, Fs, p)
+    x = x(:); noise_ref = noise_ref(:);
+    frameLen = p.frameLen; hopLen = p.hopLen; nfft = p.nfft;
+    alphaDD = p.alphaDD; Gmin = p.minGain;
+    win = hamming(frameLen, 'periodic');
+
+    [N_stft, ~, ~] = stft(noise_ref, Fs, 'Window', win, ...
+        'OverlapLength', frameLen-hopLen, 'FFTLength', nfft);
+    lambda_d = mean(abs(N_stft).^2, 2);
+    lambda_d = movmean(lambda_d, 5);
+
+    [X_stft, ~, ~] = stft(x, Fs, 'Window', win, ...
+        'OverlapLength', frameLen-hopLen, 'FFTLength', nfft);
+
+    [nFreq, nFrames] = size(X_stft);
+    Y_stft = zeros(size(X_stft));
+    G_all = zeros(nFreq, nFrames);
+    xi_prev = ones(nFreq, 1);
+
+    for k = 1:nFrames
+        X_frame = X_stft(:, k);
+        X_mag = abs(X_frame);
+        X_phase = angle(X_frame);
+        X_pow = X_mag.^2;
+        gamma_k = X_pow ./ (lambda_d + eps);
+        xi_ml = max(gamma_k - 1, 0);
+        xi = alphaDD * xi_prev + (1 - alphaDD) * xi_ml;
+        xi = max(xi, 1e-6);
+        v = xi .* gamma_k ./ (1 + xi);
+        G = zeros(nFreq,1);
+        for i = 1:nFreq
+            if v(i) < 0.01
+                G(i) = xi(i)/(1+xi(i));
+            else
+                G(i) = xi(i)/(1+xi(i)) * exp(0.5 * expint_approx(v(i)));
+            end
+        end
+        G = max(G, Gmin);
+        G = movmean(G, 3);
+        Y_mag = G .* X_mag;
+        Y_stft(:,k) = Y_mag .* exp(1j*X_phase);
+        G_all(:,k) = G;
+        xi_prev = xi;
+    end
+
+    y_out = istft(Y_stft, Fs, 'Window', win, ...
+                  'OverlapLength', frameLen-hopLen, ...
+                  'FFTLength', nfft, 'ConjugateSymmetric', true);
+    y_out = y_out(1:min(length(y_out), length(x)));
+    if length(y_out) < length(x)
+        y_out = [y_out; zeros(length(x)-length(y_out),1)];
+    end
+end
+
+
+function y = expint_approx(x)
+    if x < 1e-10
+        y = -log(1e-10) - 0.5772;
+    elseif x < 1
+        y = -log(x) - 0.5772156649;
+        term = x;
+        for n = 1:50
+            y = y + term / n;
+            term = -term * x / (n + 1);
+            if abs(term) < 1e-10, break; end
+        end
+    else
+        y = exp(-x) / x;
+        term = 1;
+        for n = 1:10
+            term = -term * n / x;
+            y = y + term * exp(-x) / x;
+            if abs(term) < 1e-10, break; end
+        end
+    end
+    y = -y;
+end

semestralka2/final.m

@@ -0,0 +1,198 @@
+clc;
+close all;
+clear;
+
+%% === LOAD AUDIO FILES ===
+[y1, Fs1] = audioread('flac.wav');
+[y2, Fs2] = audioread('flac2.wav');
+
+if size(y1,2) > 1, y1 = mean(y1,2); end
+if size(y2,2) > 1, y2 = mean(y2,2); end
+
+t1 = (0:length(y1)-1) / Fs1;
+t2 = (0:length(y2)-1) / Fs2;
+
+fprintf('Loaded flac.wav:  %d samples, Fs = %d Hz\n', length(y1), Fs1);
+fprintf('Loaded flac2.wav: %d samples, Fs = %d Hz\n', length(y2), Fs2);
+
+%% === TIME-DOMAIN WAVEFORMS ===
+figure;
+subplot(2,1,1);
+plot(t1, y1);
+xlabel('Time [s]');
+ylabel('Amplitude');
+title('flac.wav – Time Domain');
+grid on;
+
+subplot(2,1,2);
+plot(t2, y2, 'r');
+xlabel('Time [s]');
+ylabel('Amplitude');
+title('flac2.wav – Time Domain');
+grid on;
+
+%% === FREQUENCY-DOMAIN ANALYSIS (FFT) ===
+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);
+
+% --- Real / Imag components ---
+figure;
+subplot(2,1,1);
+plot(freq_shift2, real(fftshift(X2)));
+xlabel('Frequency [Hz]');
+ylabel('Real Part');
+title('flac2.wav – Real Part of FFT');
+grid on;
+subplot(2,1,2);
+plot(freq_shift2, imag(fftshift(X2)), 'r');
+xlabel('Frequency [Hz]');
+ylabel('Imag Part');
+title('flac2.wav – Imaginary Part of FFT');
+grid on;
+
+% --- Magnitude Spectrum ---
+figure;
+plot(freq_shift2, 20*log10(abs(fftshift(X2)) + eps), 'y');
+xlabel('Frequency [Hz]');
+ylabel('Magnitude [dB]');
+title('flac2.wav – Magnitude Spectrum');
+xlim([0 2000]);
+grid on;
+
+%% === PEAK DETECTION TO FIND TONAL NOISE FREQUENCIES ===
+halfN2 = floor(N2/2);
+f2 = (0:halfN2-1)*(Fs2/N2);
+mag2 = abs(X2(1:halfN2));
+
+[pks, locs] = findpeaks(mag2, 'MinPeakHeight', max(mag2)*0.1, ...
+    'MinPeakDistance', round(10*N2/Fs2));
+
+f_tonal = f2(locs);
+f_tonal = f_tonal(f_tonal > 30 & f_tonal < Fs2/2);
+f_tonal = sort(f_tonal(1:min(5, numel(f_tonal))));
+
+fprintf('\nDetected tonal noise frequencies (from flac2.wav):\n');
+fprintf('  %.1f Hz\n', f_tonal);
+
+%% === NOTCH FILTER DESIGN (Z-domain zeros/poles) ===
+r = 0.995;         % pole radius
+sos_all = [];
+g_all = 1;
+
+for i = 1:length(f_tonal)
+    f0 = f_tonal(i);
+    w0 = 2*pi*f0/Fs1;                % use Fs1 so it matches flac.wav
+    z = [exp(1j*w0) exp(-1j*w0)];    % zeros at ±ω0
+    p = r*z;                         % poles slightly inside circle
+    b = real(poly(z));
+    a = real(poly(p));
+    b = b / (sum(b)/sum(a));         % normalize gain
+    [sos_i, g_i] = tf2sos(b, a);
+    sos_all = [sos_all; sos_i];
+    g_all = g_all * g_i;
+end
+
+%% === VISUALIZE NOTCH FILTER RESPONSE ===
+figure;
+freqz(sos_all, 2048, Fs1);
+title('Notch Filter Cascade – Based on flac2 Noise Frequencies');
+
+%% === APPLY NOTCH FILTERS TO SIGNALS ===
+y1_filtered = filtfilt(sos_all, g_all, y1); % Apply to flac.wav
+y2_filtered = filtfilt(sos_all, g_all, y2); % Apply to flac2.wav
+y1_filtered = y1_filtered / max(abs(y1_filtered));
+y2_filtered = y2_filtered / max(abs(y2_filtered));
+
+audiowrite('filtered_flac.wav', y1_filtered, Fs1);
+audiowrite('filtered_flac2.wav', y2_filtered, Fs2);
+
+fprintf('\nSaved:\n');
+fprintf('  filtered_flac.wav  (filtered signal)\n');
+fprintf('  filtered_flac2.wav (filtered noise sample)\n');
+
+%% === COMPARE ORIGINAL VS FILTERED SPECTRA ===
+N_plot = length(y1);
+f_plot = (0:N_plot/2-1) * Fs1 / N_plot;
+X1_orig = fft(y1);
+X1_filt = fft(y1_filtered);
+
+figure;
+subplot(2,1,1);
+plot(f_plot, 20*log10(abs(X1_orig(1:N_plot/2))+eps));
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+title('flac.wav – Original Spectrum');
+xlim([0 2000]); grid on;
+
+subplot(2,1,2);
+plot(f_plot, 20*log10(abs(X1_filt(1:N_plot/2))+eps), 'r');
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+title('flac.wav – After Notch Filtering');
+xlim([0 2000]); grid on;
+
+%% === APPLY THE SAME FILTER TO sem2.wav ===
+fprintf('\n=== APPLYING DESIGNED FILTER TO sem2.wav ===\n');
+
+[y_sem2, Fs_sem2] = audioread('sem2.wav');
+if size(y_sem2,2) > 1, y_sem2 = mean(y_sem2,2); end
+
+if Fs_sem2 ~= Fs1
+    warning('Sample rates differ – resampling sem2.wav to match notch filter design');
+    y_sem2 = resample(y_sem2, Fs1, Fs_sem2);
+    Fs_sem2 = Fs1;
+end
+
+% --- Step 1: Apply multi‑notch filter from flac2 noise analysis ---
+y_sem2_notch = filtfilt(sos_all, g_all, y_sem2);
+
+% --- Step 2 (optional): add band‑pass filter for speech 100‑8000 Hz ---
+[z_bp, p_bp, k_bp] = butter(4, [100 8000]/(Fs_sem2/2), 'bandpass');
+[sos_bp, g_bp] = zp2sos(z_bp, p_bp, k_bp);
+y_sem2_final = filtfilt(sos_bp, g_bp, y_sem2_notch);
+
+% --- Normalize both versions ---
+y_sem2_notch = y_sem2_notch / (max(abs(y_sem2_notch)) + eps) * 0.95;
+y_sem2_final = y_sem2_final / (max(abs(y_sem2_final)) + eps) * 0.95;
+
+% --- Save results ---
+audiowrite('sem2_notch_filtered.wav', y_sem2_notch, Fs_sem2);
+audiowrite('sem2_final.wav', y_sem2_final, Fs_sem2);
+fprintf('Saved:\n');
+fprintf('  sem2_notch_filtered.wav  (notch filters only)\n');
+fprintf('  sem2_final.wav           (notch + speech band‑pass)\n');
+
+%% === VISUAL CHECK – BEFORE vs AFTER ===
+t_sem2 = (0:length(y_sem2)-1)/Fs_sem2;
+N_sem2 = length(y_sem2);
+f_sem2 = (0:N_sem2/2-1)*Fs_sem2/N_sem2;
+Y_orig = fft(y_sem2);
+Y_filt = fft(y_sem2_final);
+
+figure;
+subplot(2,1,1);
+plot(t_sem2, y_sem2);
+xlabel('Time [s]'); ylabel('Amplitude');
+title('sem2.wav – Original Time Domain'); grid on;
+
+subplot(2,1,2);
+plot(t_sem2, y_sem2_final, 'r');
+xlabel('Time [s]'); ylabel('Amplitude');
+title('sem2\_final.wav – After Notch + Band‑Pass'); grid on;
+
+figure;
+subplot(2,1,1);
+plot(f_sem2, 20*log10(abs(Y_orig(1:N_sem2/2))+eps));
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+title('sem2.wav – Original Spectrum');
+xlim([0 2000]); grid on;
+
+subplot(2,1,2);
+plot(f_sem2, 20*log10(abs(Y_filt(1:N_sem2/2))+eps), 'r');
+xlabel('Frequency [Hz]'); ylabel('Magnitude [dB]');
+title('sem2\_final.wav – After Filtering');
+xlim([0 2000]); grid on;

semestralka2/final_cleaned.m

@@ -0,0 +1,252 @@
+clc;
+close all;
+clear;
+
+% === Global figure styling (light mode defaults) ===
+set(0, 'DefaultFigureColor', 'w');        % white figure background
+set(0, 'DefaultAxesColor', 'w');          % white plotting area
+set(0, 'DefaultAxesXColor', 'k');         % black x-axis lines and text
+set(0, 'DefaultAxesYColor', 'k');         % black y-axis lines and text
+set(0, 'DefaultTextColor', 'k');          % black text for titles, labels
+set(0, 'DefaultLineColor', 'b');          % blue lines by default
+clear saveFigure
+
+% === Automatic figure export counter ===
+fig_counter = 1;
+save_fig = @saveFigure;
+
+function saveFigure(name)
+    persistent counter
+    if isempty(counter)
+        counter = 1;
+    end
+
+    % Force painters renderer and white backgrounds
+    set(gcf, 'Renderer', 'painters');
+    set(gcf, 'Color', [1 1 1]);
+    axesHandles = findall(gcf, 'Type', 'axes');
+    for ax = axesHandles'
+        set(ax, 'Color', [1 1 1]);
+    end
+
+    % Export clean EPS with no transparency problems
+    print(gcf, '-depsc2', '-r300', '-painters', sprintf('%d_%s.eps', counter, name));
+
+    counter = counter + 1;
+end
+
+[y1, Fs1] = audioread('flac.wav');
+[y2, Fs2] = audioread('flac2.wav');
+[y_sem2, Fs_sem2] = audioread('sem2.wav');
+
+% Časové priebehy
+t1 = (0:length(y1)-1)/Fs1;
+t2 = (0:length(y2)-1)/Fs2;
+t_sem2 = (0:length(y_sem2)-1)/Fs_sem2;
+
+figPos = get(0, 'DefaultFigurePosition');
+figPos(4) = figPos(4) / 2;   % polovičná výška (4 = height element)
+figure('Position', figPos);  % style defaults 
+plot(t_sem2, y_sem2);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – pôvodný (RAW) signál');
+grid on;
+save_fig('orig_sem2'); fig_counter = fig_counter + 1;
+
+figure;
+subplot(2,1,1);
+plot(t1, y1);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('flac.wav');
+grid on;
+
+subplot(2,1,2);
+plot(t2, y2, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('flac2.wav');
+grid on;
+save_fig('time'); fig_counter = fig_counter + 1;
+
+% Fourierova transformácia
+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);
+
+% Real a imag časť spektra
+figure;
+subplot(2,1,1);
+plot(freq_shift2, real(fftshift(X2)));
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac2.wav – Reálna časť spektra');
+xlim([-550 550]);
+grid on;
+
+subplot(2,1,2);
+plot(freq_shift2, imag(fftshift(X2)), 'r');
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac2.wav – Imaginárna časť spektra');
+xlim([-550 550]);
+grid on;
+save_fig('spectrum'); fig_counter = fig_counter + 1;
+
+% Real a imag časť spektra
+figure;
+subplot(2,1,1);
+plot(freq_shift1, real(fftshift(X1)));
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac.wav – Reálna časť spektra');
+xlim([-3500 3500]);
+grid on;
+
+subplot(2,1,2);
+plot(freq_shift1, imag(fftshift(X1)), 'r');
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac.wav – Imaginárna časť spektra');
+xlim([-3500 3500]);
+grid on;
+save_fig('spectrum2'); fig_counter = fig_counter + 1;
+
+% Detekcia dominantných frekvencií
+halfN2 = floor(N2 / 2);
+freq = (0:halfN2-1) * (Fs2 / N2);
+magX = abs(X2(1:halfN2));
+
+[pks, locs] = findpeaks(magX, 'NPeaks', 3, 'SortStr', 'descend');
+f_tonal = sort(freq(locs));
+
+fprintf('\nDetekovane frekvencie z (flac2.wav):\n');
+fprintf('  %.1f Hz\n', f_tonal);
+
+% SOS IIR filter
+r = 0.99;
+Fs = Fs1;
+sos = zeros(length(f_tonal), 6);
+
+for k = 1:length(f_tonal)
+    theta = 2 * pi * (f_tonal(k) / Fs);
+    b = [1, -2*cos(theta), 1];
+    a = [1, -2*r*cos(theta), r^2];
+    sos(k, :) = [b, a];
+end
+
+% Nuly a póly
+figure;
+[z, p, k_gain] = sos2zp(sos);
+zplane(z, p);
+title('Nuly a póly notch filtrov');
+xlabel('Reálna časť');
+ylabel('Imaginárna časť');
+save_fig('pz_notch'); fig_counter = fig_counter + 1;
+
+% Frekvenčná charakteristika
+figure;
+[h, w] = freqz(sos, 4096, Fs);
+freqz(sos, 4096, Fs);
+title('Frekvenčná charakteristika notch filtra');
+subplot(2,1,1); xlim([0 0.5]);
+subplot(2,1,2); xlim([0 0.5]);
+save_fig('freqz_notch'); fig_counter = fig_counter + 1;
+
+% Filtrovanie flac.wav a flac2.wav
+y1_filtered = sosfilt(sos, double(y1));
+y2_filtered = sosfilt(sos, double(y2));
+
+figure;
+subplot(2,1,1);
+plot(y2);
+title('Pôvodný signál flac2.wav');
+xlabel('Vzorky');
+ylabel('Amplitúda');
+ax1 = gca;
+grid on;
+
+subplot(2,1,2);
+plot(y2_filtered, 'r');
+title('Filtrovaný signál flac2.wav');
+xlabel('Vzorky');
+ylabel('Amplitúda');
+ax2 = gca;
+grid on;
+
+% spolocna y os pre oba grafy
+ylim_common = [min([ax1.YLim, ax2.YLim]), max([ax1.YLim, ax2.YLim])];
+ax1.YLim = ylim_common;
+ax2.YLim = ylim_common;
+save_fig('orig_vs_filtered_flac2'); fig_counter = fig_counter + 1;
+
+% Aplikácia filtrov na sem2.wav
+y_sem2_filtered = sosfilt(sos, double(y_sem2));
+
+figure;
+subplot(2,1,1);
+plot(t_sem2, y_sem2);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – pôvodný signál');
+grid on;
+
+subplot(2,1,2);
+plot(t_sem2, y_sem2_filtered, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – filtrovaný signál');
+grid on;
+
+ylim_common = [min([min(y_sem2), min(y_sem2_filtered)]), ...
+               max([max(y_sem2), max(y_sem2_filtered)])];
+subplot(2,1,1); ylim(ylim_common);
+subplot(2,1,2); ylim(ylim_common);
+save_fig('orig_vs_filtered_sem2'); fig_counter = fig_counter + 1;
+
+audiowrite('filtered_sem2.wav', y_sem2_filtered, Fs_sem2);
+
+% Butterworth filter pre rozsah ľudskej reči
+[z_bp, p_bp, k_bp] = butter(4, [100 8000]/(Fs_sem2/2), 'bandpass');
+[sos_bp, g_bp] = zp2sos(z_bp, p_bp, k_bp);
+y_sem2_final = filtfilt(sos_bp, g_bp, y_sem2_filtered);
+
+audiowrite('filtered_sem2_voiceband.wav', ...
+    y_sem2_final / max(abs(y_sem2_final)), Fs_sem2);
+
+t_sem2 = (0:length(y_sem2_filtered)-1) / Fs_sem2;
+
+figure;
+subplot(2,1,1);
+plot(t_sem2, y_sem2_filtered);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – notch filtrácia');
+grid on;
+
+subplot(2,1,2);
+plot(t_sem2, y_sem2_final, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – notch + band-pass filtrácia');
+grid on;
+
+ylim_common = [min([min(y_sem2_filtered), min(y_sem2_final)]), ...
+               max([max(y_sem2_filtered), max(y_sem2_final)])];
+subplot(2,1,1); ylim(ylim_common);
+subplot(2,1,2); ylim(ylim_common);
+save_fig('notch_vs_bandpass_sem2'); fig_counter = fig_counter + 1;
+
+% Percentuálny rozdiel band-pass filtrácie
+diff_signal = y_sem2_final - y_sem2_filtered;
+power_original = sum(y_sem2_filtered.^2);
+power_diff = sum(diff_signal.^2);
+percent_difference = (power_diff / power_original) * 100;
+fprintf('\nRozdiel medzi notch a notch+band‑pass signálom: %.2f %%\n', percent_difference);
+
+

semestralka2/final_cleaned_backup.m

@@ -0,0 +1,251 @@
+clc;
+close all;
+clear;
+
+% === Global figure styling (light mode defaults) ===
+set(0, 'DefaultFigureColor', 'w');        % white figure background
+set(0, 'DefaultAxesColor', 'w');          % white plotting area
+set(0, 'DefaultAxesXColor', 'k');         % black x-axis lines and text
+set(0, 'DefaultAxesYColor', 'k');         % black y-axis lines and text
+set(0, 'DefaultTextColor', 'k');          % black text for titles, labels
+set(0, 'DefaultLineColor', 'b');          % blue lines by default
+clear saveFigure
+
+% === Automatic figure export counter ===
+fig_counter = 1;
+save_fig = @saveFigure;
+
+function saveFigure(name)
+    persistent counter
+    if isempty(counter)
+        counter = 1;
+    end
+
+    % Force painters renderer and white backgrounds
+    set(gcf, 'Renderer', 'painters');
+    set(gcf, 'Color', [1 1 1]);
+    axesHandles = findall(gcf, 'Type', 'axes');
+    for ax = axesHandles'
+        set(ax, 'Color', [1 1 1]);
+    end
+
+    % Export clean EPS with no transparency problems
+    print(gcf, '-depsc2', '-r300', '-painters', sprintf('%d_%s.eps', counter, name));
+
+    counter = counter + 1;
+end
+
+[y1, Fs1] = audioread('flac.wav');
+[y2, Fs2] = audioread('flac2.wav');
+
+% Časové priebehy
+t1 = (0:length(y1)-1)/Fs1;
+t2 = (0:length(y2)-1)/Fs2;
+
+figure;
+subplot(2,1,1);
+plot(t1, y1);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('flac.wav');
+grid on;
+
+subplot(2,1,2);
+plot(t2, y2, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('flac2.wav');
+grid on;
+save_fig('time'); fig_counter = fig_counter + 1;
+
+% Fourierova transformácia
+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);
+
+% Real a imag časť spektra
+figure;
+subplot(2,1,1);
+plot(freq_shift2, real(fftshift(X2)));
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac2.wav – Reálna časť spektra');
+xlim([-550 550]);
+grid on;
+
+subplot(2,1,2);
+plot(freq_shift2, imag(fftshift(X2)), 'r');
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac2.wav – Imaginárna časť spektra');
+xlim([-550 550]);
+grid on;
+save_fig('spectrum'); fig_counter = fig_counter + 1;
+
+% Real a imag časť spektra
+figure;
+subplot(2,1,1);
+plot(freq_shift1, real(fftshift(X1)));
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac.wav – Reálna časť spektra');
+xlim([-3500 3500]);
+grid on;
+
+subplot(2,1,2);
+plot(freq_shift1, imag(fftshift(X1)), 'r');
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac.wav – Imaginárna časť spektra');
+xlim([-3500 3500]);
+grid on;
+save_fig('spectrum2'); fig_counter = fig_counter + 1;
+
+% Detekcia dominantných frekvencií
+halfN2 = floor(N2 / 2);
+freq = (0:halfN2-1) * (Fs2 / N2);
+magX = abs(X2(1:halfN2));
+
+[pks, locs] = findpeaks(magX, 'NPeaks', 3, 'SortStr', 'descend');
+f_tonal = sort(freq(locs));
+
+fprintf('\nDetekovane frekvencie z (flac2.wav):\n');
+fprintf('  %.1f Hz\n', f_tonal);
+
+% SOS IIR filter
+r = 0.99;
+Fs = Fs1;
+sos = zeros(length(f_tonal), 6);
+
+for k = 1:length(f_tonal)
+    theta = 2 * pi * (f_tonal(k) / Fs);
+    b = [1, -2*cos(theta), 1];
+    a = [1, -2*r*cos(theta), r^2];
+    sos(k, :) = [b, a];
+end
+
+% Nuly a póly
+figure;
+[z, p, k_gain] = sos2zp(sos);
+zplane(z, p);
+title('Nuly a póly notch filtrov');
+xlabel('Reálna časť');
+ylabel('Imaginárna časť');
+save_fig('pz_notch'); fig_counter = fig_counter + 1;
+
+% Frekvenčná charakteristika
+figure;
+[h, w] = freqz(sos, 4096, Fs);
+freqz(sos, 4096, Fs);
+title('Frekvenčná charakteristika notch filtra');
+subplot(2,1,1); xlim([0 0.5]);
+subplot(2,1,2); xlim([0 0.5]);
+save_fig('freqz_notch'); fig_counter = fig_counter + 1;
+
+% Filtrovanie flac.wav a flac2.wav
+y1_filtered = sosfilt(sos, double(y1));
+y2_filtered = sosfilt(sos, double(y2));
+
+figure;
+subplot(2,1,1);
+plot(y2);
+title('Pôvodný signál flac2.wav');
+xlabel('Vzorky');
+ylabel('Amplitúda');
+ax1 = gca;
+grid on;
+
+subplot(2,1,2);
+plot(y2_filtered, 'r');
+title('Filtrovaný signál flac2.wav');
+xlabel('Vzorky');
+ylabel('Amplitúda');
+ax2 = gca;
+grid on;
+
+ylim_common = [min([ax1.YLim, ax2.YLim]), max([ax1.YLim, ax2.YLim])];
+ax1.YLim = ylim_common;
+ax2.YLim = ylim_common;
+save_fig('orig_vs_filtered_flac2'); fig_counter = fig_counter + 1;
+
+% Aplikácia filtrov na sem2.wav
+[y_sem2, Fs_sem2] = audioread('sem2.wav');
+y_sem2_filtered = sosfilt(sos, double(y_sem2));
+
+t_sem2 = (0:length(y_sem2)-1) / Fs_sem2;
+
+figure;
+subplot(2,1,1);
+plot(t_sem2, y_sem2);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – pôvodný signál');
+grid on;
+
+subplot(2,1,2);
+plot(t_sem2, y_sem2_filtered, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – filtrovaný signál');
+grid on;
+
+ylim_common = [min([min(y_sem2), min(y_sem2_filtered)]), ...
+               max([max(y_sem2), max(y_sem2_filtered)])];
+subplot(2,1,1); ylim(ylim_common);
+subplot(2,1,2); ylim(ylim_common);
+save_fig('orig_vs_filtered_sem2'); fig_counter = fig_counter + 1;
+
+audiowrite('filtered_sem2.wav', y_sem2_filtered, Fs_sem2);
+
+% Butterworth filter pre rozsah ľudskej reči
+[z_bp, p_bp, k_bp] = butter(4, [100 8000]/(Fs_sem2/2), 'bandpass');
+[sos_bp, g_bp] = zp2sos(z_bp, p_bp, k_bp);
+y_sem2_final = filtfilt(sos_bp, g_bp, y_sem2_filtered);
+
+audiowrite('filtered_sem2_voiceband.wav', ...
+    y_sem2_final / max(abs(y_sem2_final)), Fs_sem2);
+
+t_sem2 = (0:length(y_sem2_filtered)-1) / Fs_sem2;
+
+figure;
+subplot(2,1,1);
+plot(t_sem2, y_sem2_filtered);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – notch filtrácia');
+grid on;
+
+subplot(2,1,2);
+plot(t_sem2, y_sem2_final, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – notch + band-pass filtrácia');
+grid on;
+
+ylim_common = [min([min(y_sem2_filtered), min(y_sem2_final)]), ...
+               max([max(y_sem2_filtered), max(y_sem2_final)])];
+subplot(2,1,1); ylim(ylim_common);
+subplot(2,1,2); ylim(ylim_common);
+save_fig('notch_vs_bandpass_sem2'); fig_counter = fig_counter + 1;
+
+% Percentuálny rozdiel band-pass filtrácie
+diff_signal = y_sem2_final - y_sem2_filtered;
+power_original = sum(y_sem2_filtered.^2);
+power_diff = sum(diff_signal.^2);
+percent_difference = (power_diff / power_original) * 100;
+fprintf('\nRozdiel medzi notch a notch+band‑pass signálom: %.2f %%\n', percent_difference);
+
+figPos = get(0, 'DefaultFigurePosition');
+figPos(4) = figPos(4) / 2;   % polovičná výška (4 = height element)
+figure('Position', figPos);  % použije globálne defaulty štýlu
+plot(t_sem2, y_sem2);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('sem2.wav – pôvodný (RAW) signál');
+grid on;
+save_fig('orig_sem2'); fig_counter = fig_counter + 1;
+

semestralka2/spd_sumu.m

@@ -0,0 +1,207 @@
+clc;
+close all;
+clear;
+
+% === Global figure styling (light mode defaults) ===
+set(0, 'DefaultFigureColor', 'w');        % white figure background
+set(0, 'DefaultAxesColor', 'w');          % white plotting area
+set(0, 'DefaultAxesXColor', 'k');         % black x-axis lines and text
+set(0, 'DefaultAxesYColor', 'k');         % black y-axis lines and text
+set(0, 'DefaultTextColor', 'k');          % black text for titles, labels
+set(0, 'DefaultLineColor', 'b');          % blue lines by default
+clear saveFigure
+
+% === Automatic figure export counter ===
+fig_counter = 1;
+save_fig = @saveFigure;
+
+function saveFigure(name)
+    persistent counter
+    if isempty(counter)
+        counter = 1;
+    end
+
+    % Force painters renderer and white backgrounds
+    set(gcf, 'Renderer', 'painters');
+    set(gcf, 'Color', [1 1 1]);
+    axesHandles = findall(gcf, 'Type', 'axes');
+    for ax = axesHandles'
+        set(ax, 'Color', [1 1 1]);
+    end
+
+    % Export clean EPS with no transparency problems
+    print(gcf, '-depsc2', '-r300', '-painters', sprintf('%d_%s.eps', counter, name));
+
+    counter = counter + 1;
+end
+
+[y1, Fs1] = audioread('flac.wav');
+[y2, Fs2] = audioread('flac2.wav');
+
+% Časové priebehy
+t1 = (0:length(y1)-1)/Fs1;
+t2 = (0:length(y2)-1)/Fs2;
+
+figure;
+subplot(2,1,1);
+plot(t1, y1);
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('flac.wav');
+grid on;
+
+subplot(2,1,2);
+plot(t2, y2, 'r');
+xlabel('Čas [s]');
+ylabel('Amplitúda');
+title('flac2.wav');
+grid on;
+save_fig('time'); fig_counter = fig_counter + 1;
+
+% Fourierova transformácia
+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);
+
+% Real a imag časť spektra
+figure;
+subplot(2,1,1);
+plot(freq_shift2, real(fftshift(X2)));
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac2.wav – Reálna časť spektra');
+xlim([-550 550]);
+grid on;
+
+subplot(2,1,2);
+plot(freq_shift2, imag(fftshift(X2)), 'r');
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac2.wav – Imaginárna časť spektra');
+xlim([-550 550]);
+grid on;
+save_fig('spectrum'); fig_counter = fig_counter + 1;
+
+% Real a imag časť spektra
+figure;
+subplot(2,1,1);
+plot(freq_shift1, real(fftshift(X1)));
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac.wav – Reálna časť spektra');
+xlim([-3500 3500]);
+grid on;
+
+subplot(2,1,2);
+plot(freq_shift1, imag(fftshift(X1)), 'r');
+xlabel('Frekvencia [Hz]');
+ylabel('Amplitúda');
+title('flac.wav – Imaginárna časť spektra');
+xlim([-3500 3500]);
+grid on;
+save_fig('spectrum2'); fig_counter = fig_counter + 1;
+
+% Určenie typu šumu na základe spektrálneho sklonu
+[PSD, f] = pwelch(y1, hamming(4096), [], [], Fs1, 'onesided'); % PSD
+mask = f > 10 & f < Fs1/2;
+logf = log10(f(mask));
+logP = log10(PSD(mask));
+p = polyfit(logf, logP, 1);
+slope = p(1);
+intercept = p(2);
+fprintf('Spektrálny sklon (alfa) = %.2f\n', slope);
+
+% Vizualizácia
+figure;
+loglog(f, PSD); hold on;
+loglog(f(mask), 10.^(polyval(p, logf)), 'r', 'LineWidth', 1.5);
+xlabel('Frekvencia [Hz]');
+ylabel('Výkonová spektrálna hustota (PSD)');
+title(sprintf('Spektrálny sklon = %.2f', slope));
+lgd = legend('Meraná PSD', 'Lineárna aproximácia (log-log)', 'Location', 'best');
+set(lgd, 'Color', 'w', 'TextColor', 'k');
+grid on;
+save_fig('spectral_slope'); fig_counter = fig_counter + 1;
+
+% === IIR Notch filtering flac2.wav and noise type analysis ===
+
+fprintf('\n=== IIR Notch Filtering of flac2.wav ===\n');
+
+% Detekcia dominantných frekvencií v flac2.wav pre notch filtráciu
+halfN = floor(N2 / 2);
+freq_pos = (0:halfN - 1) * (Fs2 / N2);
+magX2 = abs(X2(1:halfN));
+[pks2, locs2] = findpeaks(magX2, 'NPeaks', 3, 'SortStr', 'descend');
+f_notch = sort(freq_pos(locs2));
+
+fprintf('\nDetekované dominantné frekvencie vo flac2.wav:\n');
+fprintf('  %.1f Hz\n', f_notch);
+
+% === Konštrukcia a aplikácia IIR notch filtra ===
+r = 0.98;  % ostrosť zárezu (0.98 až 0.995)
+Fs = Fs2;
+sos_notch = zeros(length(f_notch), 6);
+
+for k = 1:length(f_notch)
+    theta = 2 * pi * (f_notch(k) / Fs);
+    b = [1, -2*cos(theta), 1];
+    a = [1, -2*r*cos(theta), r^2];
+    sos_notch(k, :) = [b, a];
+end
+
+% Aplikácia filtra na flac2.wav
+y2_notched = sosfilt(sos_notch, double(y2));
+
+% Uloženie výsledného signálu
+audiowrite('flac2_notched.wav', y2_notched / max(abs(y2_notched)), Fs2);
+
+% Zobrazenie pred/po filtrácii
+figure;
+subplot(2,1,1);
+plot(y2);
+xlabel('Vzorky'); ylabel('Amplitúda');
+title('Pôvodný signál flac2.wav');
+grid on;
+
+subplot(2,1,2);
+plot(y2_notched, 'r');
+xlabel('Vzorky'); ylabel('Amplitúda');
+title('Filtrovaný signál – flac2.wav (IIR Notch)');
+grid on;
+save_fig('flac2_notched'); fig_counter = fig_counter + 1;
+
+% === Analýza zvyškového šumu po IIR filtri ===
+[PSD_res, f_res] = pwelch(y2_notched, hamming(4096), [], [], Fs2, 'onesided');
+mask_res = f_res > 10 & f_res < Fs2/2;
+logf_res = log10(f_res(mask_res));
+logP_res = log10(PSD_res(mask_res));
+p_res = polyfit(logf_res, logP_res, 1);
+slope_res = p_res(1);
+
+fprintf('\nSpektrálny sklon po IIR filtrácii (alfa) = %.2f\n', slope_res);
+
+% Vizualizácia
+figure;
+loglog(f_res, PSD_res); hold on;
+loglog(f_res(mask_res), 10.^(polyval(p_res, logf_res)), 'r', 'LineWidth', 1.5);
+xlabel('Frekvencia [Hz]');
+ylabel('Výkonová spektrálna hustota (PSD)');
+title(sprintf('Zvyškový šum po IIR notch filtrácii (α = %.2f)', slope_res));
+legend('PSD po IIR filtri', 'Lineárna aproximácia', 'Location', 'best');
+grid on;
+save_fig('residual_noise_flac2'); fig_counter = fig_counter + 1;
+
+% Interpretácia typu šumu podľa alfa
+if slope_res > -0.5
+    noise_type = 'biely šum';
+elseif slope_res > -1.5
+    noise_type = 'ružový šum';
+else
+    noise_type = 'hnedý šum';
+end
+
+fprintf('Typ zostávajúceho šumu po notch filtrácii: %s\n', noise_type);