commit eb210b9f462b279f5ae0887d23d52c0c4b7fced0
parent 2d0e20368adafedb0ce9a7df5ad56e1fda5e61d5
Author: Martin Kloeckner <mjkloeckner@gmail.com>
Date: Thu, 27 Nov 2025 16:05:18 -0300
updated `tp/*`
added `tp/plot/potencia_db_espectro_ventana_comparacion_hamming_rect.png` and
modified scripts accordingly to generate it
Diffstat:
5 files changed, 220 insertions(+), 38 deletions(-)
diff --git a/tp/plot/potencia_db_espectro_ventana_comparacion_hamming_rect.png b/tp/plot/potencia_db_espectro_ventana_comparacion_hamming_rect.png
Binary files differ.
diff --git a/tp/scripts/data.py b/tp/scripts/data.py
@@ -36,10 +36,14 @@ file2_filter2_output = np.convolve(file2_data, filter2_h, mode='same')
# todas las canciones tienen formato mp3 y 44100Hz de frecuencia de muestreo
canciones_dataset_dir = data_dir + 'canciones/'
canciones_dataset = []
+canciones_dataset_common_fs = 44100
-for filename in os.listdir(canciones_dataset_dir):
+for i, filename in enumerate(sorted(os.listdir(canciones_dataset_dir))):
filename = canciones_dataset_dir + filename
print(filename)
file_data, file_fs = librosa.load(data_dir + filename, sr=None, mono=True)
canciones_dataset.append(file_data)
file_fs = int(file_fs)
+
+ # cargar solo una canción
+ break
diff --git a/tp/scripts/main.py b/tp/scripts/main.py
@@ -34,7 +34,16 @@ def segunda_parte():
a4_violin_cutoff()
def tercera_parte():
- filtro_fir_deducido()
+ # h, fs = filtro_fir_deducido()
+
+ # filtro_fir_analisis(h, fs)
+
+ # for i in len(canciones_dataset):
+ # filtro_fir_filtrar_comparar_espectogramas(
+ # h, canciones_dataset[i], canciones_dataset_common_fs)
+
+ # analisis_freq_ventanas()
+
# primera_parte()
# segunda_parte()
diff --git a/tp/scripts/tercera_parte.py b/tp/scripts/tercera_parte.py
@@ -1,5 +1,6 @@
from utils import *
from scipy.signal import firwin, freqz, tf2zpk, get_window
+from scipy.fft import fftshift
############################# Tercera parte ###################################
@@ -7,6 +8,8 @@ cutoff = 2650 # frecuencia de corte
M = 700 # orden FIR (número de coeficientes)
def filtro_fir():
+ fs = 44100
+
# Diseño FIR pasabajos con ventana
b = firwin(M, cutoff, fs=fs, window='hamming')
@@ -51,12 +54,12 @@ def filtro_fir_polos_y_ceros(a):
save_plot(fig, "polos_y_ceros_pasa-bajos_fir")
-fs = 44100
-
# a diferencia de la funcion `filtro_fir` se genera el filtro mediante
# operaciones elementales, como la multiplicacion por ventana, en lugar de usar
# una funcion de libreria externa como `firwin`
def filtro_fir_deducido():
+ fs = 44100
+
# respuesta ideal pasabajos: sinc centrada en M/2
n = np.arange(M + 1)
wc = 2*np.pi*cutoff / fs
@@ -65,17 +68,15 @@ def filtro_fir_deducido():
# se normaliza la ganancia a 1 multiplicando por 2.0*(fc/fs)
h_ideal = np.sinc(2.0 * (cutoff/fs) * (n - M/2))
- # ventana de Hamming
+ # ventana de hamming (de acuerdo a la formula de wikipedia)
v_hamming = 0.54 - 0.46 * np.cos(2*np.pi*n/M)
- v_rectangular = [1 if i < 350 else 0 for i in range(0, M)]
-
# respuesta del filtro FIR (version acotada de la sinc)
h = h_ideal * v_hamming
# se normaliza para tener ganancia unitaria para frecuancias <= fc
h = h / np.sum(h)
- return h
+ return h, fs
def filtro_fir_analisis(h, fs):
fig, ax = dtime_plot(M, h, "respuesta_al_impulso_filtro_fir",
@@ -91,34 +92,114 @@ def filtro_fir_analisis(h, fs):
# polos y ceros
filtro_fir_polos_y_ceros(h)
+def filtro_fir_filtrar_comparar_espectogramas(h, data, fs):
+ spectogram_plot(fs, data,
+ f"filtro_fir_espectograma_original_44100Hz",
+ N=1024, ylim=[0, 20000])
+
+ filter_output = np.convolve(data, h, mode='same')
+ spectogram_plot(fs, filter_output,
+ f"filtro_fir_espectograma_filtrada_{cutoff}Hz",
+ t=0, N=1024, ylim=[0, 20000])
+
+
+def analisis_freq_ventanas():
+ N = 1024
+
+ v_rectangular = np.ones(N)
+ v_hamming = np.hamming(N)
+
+ # freq_plot(44100, v_hamming, "v_hamming_freq", f_max=2000, N=N*8)
+ # freq_plot(44100, v_rectangular, "v_rectangular_freq", f_max=2000, N=N*8)
+
+ v_rect_fft, _ = freq_compute_fft(44100, v_rectangular, N=N*8)
+ v_hamm_fft, _ = freq_compute_fft(44100, v_hamming, N=N*8)
+
+ # se centra el lobulo principal (frecuencia 0 en el centro del arreglo)
+ v_rect_fft = fftshift(v_rect_fft)
+ v_hamm_fft = fftshift(v_hamm_fft)
+
+ v_rect_potencia = np.abs(v_rect_fft)**2
+ v_hamm_potencia = np.abs(v_hamm_fft)**2
+
+ # se agrega 1e-12 para evitar dividir por cero
+ v_rect_potencia_db = 10 * np.log10(v_rect_potencia + 1e-12)
+ v_hamm_potencia_db = 10 * np.log10(v_hamm_potencia + 1e-12)
+
+ # Normalizar para que el pico del lóbulo principal sea 0 dB
+ v_rect_potencia_db = v_rect_potencia_db - np.max(v_rect_potencia_db)
+ v_hamm_potencia_db = v_hamm_potencia_db - np.max(v_hamm_potencia_db)
+
+ # Gráfico de la Ventana Rectangular
+ f = np.linspace(-0.5, 0.5, N*8)
+
+ fig, ax = freq_graph_data_norm(f, v_rect_potencia_db,
+ x_min=-0.04, x_max=0.04, y_min=-200, y_max=10,
+ show=False)
+ save_plot(fig, "potencia_db_espectro_ventana_rectangular")
+
+ fig, ax = freq_graph_data_norm(f, v_hamm_potencia_db,
+ x_min=-0.04, x_max=0.04, y_min=-200, y_max=10,
+ show=False)
+ save_plot(fig, "potencia_db_espectro_ventana_hamming")
+
+
+ data_arr = [v_rect_potencia_db, v_hamm_potencia_db]
+ leg_arr = ["Potencia espectro ventana rectangular",
+ "Potencia espectro ventana Hamming"]
+
+ fig, ax = freq_graph_multiple_data_norm(f, data_arr, leg_arr,
+ x_min=-0.04, x_max=0.04, y_min=-100, y_max=40,
+ show=False)
+
+ save_plot(fig, "potencia_db_espectro_ventana_comparacion_hamming_rect")
+
+ # plt.figure(figsize=(10, 6))
+ # plt.plot(f, v_rect_potencia_db,
+ # label='Ventana Rectangular', linewidth=2, color='tab:blue')
+
+ # Gráfico de la Ventana de Hamming
+ # plt.plot(f, v_hamm_potencia_db,
+ # label='Ventana de Hamming', linewidth=2, color='tab:orange')
+
+ # Enfocarse en el lóbulo principal
+ # plt.xlim(-0.1, 0.1)
+
+ # Configuración del Gráfico
+ # plt.title('Comparación de Espectros de Potencia de Ventanas (Escala Logarítmica)')
+ # plt.xlabel('Frecuencia Normalizada (f/Fs)')
+ # plt.ylabel('Potencia Normalizada (dB)')
+
+ # Limitar el eje Y para ver bien los lóbulos laterales
+ # plt.ylim(-100, 5)
+ # Limitar el eje X para enfocarse en el lóbulo principal
+ # plt.xlim(-0.1, 0.1)
-"""
-def ejemplo_cancion_filtrado_con_filtro_fir():
- # spectogram_plot(file_fs, file_data,
- # f"espectograma_fs_original_44100Hz", N=1024, ylim=[0, 20000])
+ # plt.grid(True, which="both", linestyle='--', alpha=0.7)
+ # plt.legend()
+ # plt.show()
- file_filter_output = np.convolve(file_data, h, mode='same')
+ # freq_plot(44100, v_hamming_potencia, "v_hamming_potencia", f_max=20000, N=8192)
+ # freq_plot(44100, v_rectangular_potencia, "v_rectangular_potencia", f_max=20000, N=8192)
- # spectogram_plot(file_fs, file_filter_output,
- # f"espectograma_fs_{cutoff}Hz", t=0, N=1024, ylim=[0, 20000])
+ # interval_fft, interval_freqs = freq_compute_fft(fs, data, t, dt, N=N)
+ # interval_fft, interval_freqs = freq_compute_fft(fs, data, t, dt, N=N)
- freq_plot(44100, v_hamming, "v_hamming_freq", f_max=8000)
- freq_plot(44100, v_rectangular, "v_rectangular_freq", f_max=8000)
+ # f = np.linspace(-0.5, 0.5, n_fft)
- freq_compute_fft(44100, v_hamming)
+ # filter_output = np.convolve(data, h, mode='same')
- for i in [512, 1024, 2048]:
- # for window in ['boxcar', 'bartlett', 'hamming']:
+ # for i in [512, 1024, 2048]:
+ # # for window in ['boxcar', 'bartlett', 'hamming']:
- # spectogram_plot(file1_fs, file_filter_output,
- # f"espectograma_submuestreado_{window}_{i:04d}", N=i,
- # win=window, ylim=[0, 3000], t=5, dt=1)
+ # # spectogram_plot(file1_fs, file_filter_output,
+ # # f"espectograma_submuestreado_{window}_{i:04d}", N=i,
+ # # win=window, ylim=[0, 3000], t=5, dt=1)
- N = 1024
- beta = 8.6
- kaiser_window = get_window(("kaiser", beta), N)
- spectogram_plot(file1_fs, file_filter_output,
- f"espectograma_submuestreado_kaiser_window_{i:04d}", N=i,
- win=kaiser_window, ylim=[0, 3000], t=5, dt=1)
-"""
+ # N = 1024
+ # beta = 8.6
+ # kaiser_window = get_window(("kaiser", beta), N)
+ # spectogram_plot(canciones_dataset_common_fs, filter_output,
+ # f"espectograma_submuestreado_kaiser_window_{i:04d}", N=i,
+ # win=kaiser_window, ylim=[0, 3000], t=5, dt=1)
diff --git a/tp/scripts/utils.py b/tp/scripts/utils.py
@@ -128,17 +128,28 @@ def time_plot(fs, data, save_name="", t=0, dt=0, a=0, da=0):
return fig, ax
-def save_plot(fig, name):
+def save_plot(fig, name, overwrite=True):
base_name = os.path.basename(name)
file_name, ext = os.path.splitext(base_name)
file_path_no_ext = f'{plot_dir}{file_name}'
save_name = f'{file_path_no_ext}.png'
+
+ if overwrite == False:
+ if os.path.exists(save_name):
+ i = 1
+ while True:
+ new_save_name = f'{file_path_no_ext}_{i:02d}'
+ if not os.path.exists(f'{new_save_name}.png'):
+ save_name = f'{new_save_name}.png'
+ break
+ i += 1
+
print(save_name)
# crea carpeta para plots
os.makedirs(plot_dir, exist_ok=True)
- fig.savefig(save_name, dpi=250, bbox_inches="tight")
+ fig.savefig(save_name, dpi=100, bbox_inches="tight")
plt.close(fig) # liberar memoria
def save_to_wav(fs, data, save_name):
@@ -168,7 +179,7 @@ def freq_graph_multiple_data(data, f_min=0, f_max=0, y_min=0, y_max=0, show=True
y = np.abs(fft[:N // 2])
axis.plot(x, y, label=label, alpha=0.90, linewidth=((len(data)-i-1)*0.5 + 1.5))
- axis.set(xlabel='Frecuencia [Hz]', ylabel='Magnitud')
+ axis.set(xlabel='Frecuencia [Hz]', ylabel='Magnitud [dB]')
axis.minorticks_on()
axis.grid(True, which='major', color='black', linestyle=':', linewidth=1.00)
@@ -197,12 +208,89 @@ def freq_graph_multiple_data(data, f_min=0, f_max=0, y_min=0, y_max=0, show=True
return fig, axis
+def freq_graph_data_norm(x, y, x_min=0, x_max=0, y_min=0, y_max=0, show=True,
+ xlabel="", ylabel=""):
+
+ fig, axis = plt.subplots(figsize=(8, 4))
+
+ axis.plot(x, y)
+
+ if xlabel == "":
+ axis.set(xlabel='Frecuencia normalizada')
+
+ if ylabel == "":
+ axis.set(ylabel='Magnitud [dB]')
+
+ axis.minorticks_on()
+ axis.grid(True, which='major', color='black', linestyle=':', linewidth=1.00)
+ axis.grid(True, which='minor', color='black', linestyle=':', linewidth=0.50)
+
+ # configuracion de ticks del eje x
+ axis.xaxis.set_major_locator(MaxNLocator(nbins=5))
+ axis.xaxis.set_minor_locator(AutoMinorLocator(5))
+
+ axis.yaxis.set_major_locator(MaxNLocator(nbins=5))
+ axis.yaxis.set_minor_locator(AutoMinorLocator(4))
+
+ plt.tight_layout()
+
+ axis.set_xlim([x_min, x_max if x_max != 0 else 1])
+ axis.set_ylim([y_min, y_max if y_max != 0 else 1.05*max(y)])
+
+ if show:
+ plt.show()
+
+ return fig, axis
+
+def freq_graph_multiple_data_norm(x, y_arr, labels,
+ x_min=0, x_max=0, y_min=0, y_max=0,
+ show=True, xlabel="", ylabel=""):
+
+ fig, axis = plt.subplots(figsize=(8, 4))
+
+ for i, y in enumerate(y_arr):
+ axis.plot(x, y, label=labels[i], alpha=0.75)
+
+ if xlabel == "":
+ axis.set(xlabel='Frecuencia normalizada')
+
+ if ylabel == "":
+ axis.set(ylabel='Magnitud [dB]')
+
+ axis.minorticks_on()
+ axis.grid(True, which='major', color='black', linestyle=':', linewidth=1.00)
+ axis.grid(True, which='minor', color='black', linestyle=':', linewidth=0.50)
+
+ # configuracion de ticks del eje x
+ axis.xaxis.set_major_locator(MaxNLocator(nbins=5))
+ axis.xaxis.set_minor_locator(AutoMinorLocator(5))
+
+ axis.yaxis.set_major_locator(MaxNLocator(nbins=5))
+ axis.yaxis.set_minor_locator(AutoMinorLocator(4))
+
+ plt.tight_layout()
+ axis.legend(loc='upper right')
+
+ axis.set_xlim([x_min, x_max if x_max != 0 else 1])
+ axis.set_ylim([y_min, y_max if y_max != 0 else 1.05*max(y)])
+
+ if show:
+ plt.show()
+
+ return fig, axis
+
-def freq_graph_data(x, y, f_min=0, f_max=0, y_min=0, y_max=0, show=True):
+def freq_graph_data(x, y, f_min=0, f_max=0, y_min=0, y_max=0, show=True,
+ xlabel="", ylabel=""):
fig, axis = plt.subplots(figsize=(8, 4))
axis.plot(x, y)
- axis.set(xlabel='Frecuencia [Hz]', ylabel='Magnitud')
+
+ if xlabel == "":
+ axis.set(xlabel='Frecuencia [Hz]')
+
+ if ylabel == "":
+ axis.set(ylabel='Magnitud')
axis.minorticks_on()
axis.grid(True, which='major', color='black', linestyle=':', linewidth=1.00)
@@ -249,9 +337,9 @@ def freq_compute_fft(fs, data, t=0, dt=0, N=0):
# hace la transformacion a frecuencias y pasa lo transformado a `freq_graph_data`
def freq_plot(fs, data, save_name="", f_min=0, f_max=0, y_min=0, y_max=0,
- t=0, dt=0, a=0, da=0, show=False):
+ t=0, dt=0, a=0, da=0, show=False, N=0):
- interval_fft, interval_freqs = freq_compute_fft(fs, data, t, dt)
+ interval_fft, interval_freqs = freq_compute_fft(fs, data, t, dt, N=N)
N = len(interval_fft)
# se toma la parte positiva en ambos casos (primer parte del arreglo)
@@ -312,7 +400,7 @@ def spectogram_plot(fs, data, save_name="", t=0, dt=0, N=1024, overlp=16, win='h
plt.show()
else:
if save_name != "":
- save_plot(fig, save_name)
+ save_plot(fig, save_name, overwrite=False)
return fig, axis
