Adquisición de audio y FFT en MATLAB

Actualizado:

Capturar audio desde el micrófono y ver sus frecuencias con la FFT (transformada rápida de Fourier) es una práctica clásica en los cursos de señales. Con ella puedes descubrir la frecuencia de un silbido, la nota de una guitarra o el zumbido de 60 Hz de la red eléctrica. En este post te explico el código original que usé en 2012, por qué hoy ya no funciona en versiones nuevas de MATLAB, y una versión moderna y corregida, con la explicación de cada paso y un ejemplo numérico.

El código original (histórico)

Este programa usaba la interfaz antigua de la Data Acquisition Toolbox (analoginput) para leer la tarjeta de sonido de Windows ('winsound'), graficar la señal mientras se grababa y luego calcular su espectro:

clc; close all; clear all;
AI = analoginput('winsound');
addchannel(AI,1);
duration = 10;                          % adquisición de 10 segundos
ActualRate = get(AI,'SampleRate');
set(AI,'SamplesPerTrigger',duration*ActualRate)
figure
P = plot(zeros(1000,1));
T = title([sprintf('Peekdata calls: '), num2str(0)]);
xlabel('Samples'), axis([0 60000 -1 1]), grid on
start(AI)
i = 1;
while AI.SamplesAcquired < AI.SamplesPerTrigger
    while AI.SamplesAcquired < 1000*i
    end
    data = peekdata(AI,60000);
    set(P,'ydata',data);
    set(T,'String',[sprintf('Peekdata calls: '),num2str(i)]);
    drawnow
    i = i + 1;
end
wait(AI,2)
delete(AI)
clear AI
T = 1/8000; N = 60000; fs = 1/T;
H_data = abs(fft(data));                % magnitud del espectro
freq = (1:N/2)*fs/N;
plot(freq,H_data(1:N/2),'k');
xlim([0 100])
%% by Luis Carlos Portilla

Esa interfaz "legacy" (analoginput, addchannel, peekdata) fue retirada en las versiones modernas de MATLAB, así que el código produce error. Además tenía tres detalles a corregir:

  • Asumía fs = 8000 Hz, aunque la frecuencia real era la que devolvía ActualRate.
  • El eje de frecuencias empezaba en 1·fs/N en lugar de 0, lo que desplazaba todo el espectro un punto.
  • xlim([0 100]) solo mostraba hasta 100 Hz, por lo que la voz o un silbido (cientos o miles de Hz) no se veían.

Versión moderna y corregida

clc; clear; close all;

fs = 8000;          % frecuencia de muestreo (Hz)
dur = 5;            % segundos de grabación
rec = audiorecorder(fs, 16, 1);

disp('Grabando... silba o habla cerca del micrófono');
recordblocking(rec, dur);
x = getaudiodata(rec);          % vector columna entre -1 y 1
N = length(x);

% --- Señal en el tiempo ---
t = (0:N-1)/fs;
subplot(2,1,1); plot(t, x); grid on;
xlabel('Tiempo (s)'); ylabel('Amplitud'); title('Señal capturada');

% --- Espectro de magnitud (un solo lado) ---
X = abs(fft(x)) / N;
X = X(1:floor(N/2)+1);
X(2:end-1) = 2*X(2:end-1);
f = (0:floor(N/2)) * fs / N;     % eje de frecuencias desde 0 Hz

subplot(2,1,2); plot(f, X, 'k'); grid on;
xlabel('Frecuencia (Hz)'); ylabel('|X(f)|'); title('Espectro (FFT)');
xlim([0 fs/2]);

[~, k] = max(X(2:end));          % ignorar la componente de 0 Hz
fprintf('Frecuencia dominante: %.1f Hz\n', f(k+1));

Qué hace cada parte

  1. audiorecorder(fs, 16, 1) prepara una grabación a 8000 Hz, 16 bits, mono. Funciona con el micrófono del computador sin toolboxes adicionales.
  2. recordblocking graba durante dur segundos y espera a que termine.
  3. fft(x) calcula la transformada. Dividimos entre N para obtener amplitudes en la misma escala que la señal.
  4. Nos quedamos con la primera mitad del espectro (la otra mitad es simétrica para señales reales) y duplicamos los valores intermedios para conservar la amplitud.
  5. f = (0:N/2)·fs/N construye el eje de frecuencias desde 0 Hz hasta fs/2 (frecuencia de Nyquist).
  6. max encuentra el pico más alto: la frecuencia dominante.

Ejemplo numérico: resolución en frecuencia

La separación entre dos puntos del espectro es Δf = fs / N.

  • Con fs = 8000 Hz y 5 s de grabación: N = 40 000 muestras → Δf = 8000 / 40 000 = 0,2 Hz.
  • En el código original: N = 60 000 a 8000 Hz → Δf ≈ 0,13 Hz.
  • La frecuencia máxima visible es fs/2 = 4000 Hz.

Grabar más tiempo mejora la resolución (puntos más juntos); aumentar fs permite ver frecuencias más altas.

Prueba sin micrófono

Para comprobar que el código funciona, reemplaza la grabación por una señal conocida:

fs = 8000; t = (0:1/fs:2-1/fs)';
x = 0.8*sin(2*pi*440*t) + 0.3*sin(2*pi*1000*t) + 0.05*randn(size(t));

El espectro debe mostrar dos picos, en 440 Hz (amplitud ≈ 0,8) y en 1000 Hz (≈ 0,3), sobre un piso de ruido bajo. Si quieres generar otras señales de prueba, mira onda senoidal en MATLAB.

Errores comunes

  • "Undefined function 'analoginput'": estás usando el código antiguo en una versión moderna. Usa audiorecorder.
  • El espectro sale plano o vacío: el micrófono no está grabando. Revisa los permisos de micrófono de Windows y prueba con sound(x, fs).
  • Aparece un pico enorme en 0 Hz: es la componente continua. Réstale la media a la señal: x = x - mean(x).
  • Picos en 60 Hz y sus múltiplos: es interferencia de la red eléctrica (en Colombia, 60 Hz). Aleja el micrófono de cables y transformadores.

Preguntas frecuentes

¿Por qué solo se grafica la mitad del espectro?

Para señales reales, la segunda mitad de la FFT es un reflejo de la primera y no aporta información nueva.

¿Qué frecuencia de muestreo uso para música?

44 100 Hz o 48 000 Hz, que permiten ver hasta 22 o 24 kHz. Lo explico en resolución de bits y muestreo de audio en MATLAB.

¿Qué es mejor, fft o pwelch?

fft es ideal para aprender. Para señales ruidosas, pwelch (Signal Processing Toolbox) promedia varios segmentos y da un espectro más limpio.

Conclusión

Con audiorecorder y fft puedes capturar audio y ver sus frecuencias en pocas líneas. Recuerda construir el eje de frecuencias desde 0 Hz con fs/N, mostrar hasta fs/2 y elegir la duración según la resolución que necesites.

Comentarios

No hay comentarios

Suscribirse a: Enviar comentarios (Atom)