audio - vibracion - ¿Cómo encontrar la frecuencia fundamental de un sonido de cuerda de guitarra?
ondas estacionarias (3)
Quiero construir una aplicación de afinador de guitarra para Iphone. Mi objetivo es encontrar la frecuencia fundamental del sonido generado por una cuerda de guitarra. He utilizado bits de código de la muestra aurioTouch proporcionada por Apple para calcular el espectro de frecuencias y encuentro la frecuencia con la amplitud más alta. Funciona bien para sonidos puros (los que solo tienen una frecuencia) pero para sonidos de una cuerda de guitarra produce resultados incorrectos. He leído que esto se debe a los armónicos generados por la cuerda de guitarra que pueden tener amplitudes más altas que la fundamental. ¿Cómo puedo encontrar la frecuencia fundamental para que funcione con las cuerdas de guitarra? ¿Hay una biblioteca de código abierto en C / C ++ / Obj-C para el análisis de sonido (o procesamiento de señales)?
Encontrar los tonos musicales en un acorde es mucho más difícil que estimar el tono de una sola cuerda o nota tocada a la vez. Los armónicos de las múltiples notas en un acorde pueden estar superpuestos y entrelazados. Y todas las notas en acordes comunes pueden estar en frecuencias de armónicos para una o más notas de tono más bajo que no existen.
Para notas simples, la autocorrelación es una técnica común utilizada por algunos afinadores de guitarra. Pero con la autocorrelación, debe tener en cuenta una posible incertidumbre de octava, ya que las guitarras pueden producir connotaciones inarmónicas y en descomposición que, por lo tanto, no coinciden exactamente entre el período de tono y el período de tono. Cepstrum y Harmonic Product Spectrum son otros dos métodos de estimación de tono que pueden o no tener diferentes problemas, dependiendo de la guitarra y la nota.
RAPT parece ser un algoritmo publicado para una estimación de tono más robusta. YIN es otro.
Además, Objective C es un superconjunto de ANSI C. Por lo tanto, puede usar cualquier rutina C DSP que encuentre para la estimación de tono dentro de una aplicación de Objective C.
Puede utilizar la autocorrelación de la señal, que es la transformada inversa de la magnitud al cuadrado de la DFT. Si está muestreando a 44100 muestras / s, entonces un fundamental de 82.4 Hz es aproximadamente 535 muestras, mientras que 1479.98 Hz es aproximadamente 30 muestras. Busque el retraso positivo máximo en ese rango (por ejemplo, de 28 a 560). Asegúrese de que su ventana tenga al menos dos períodos del fundamental más largo, que serían 1070 muestras aquí. Para la siguiente potencia de dos es un búfer de 2048 muestras. Para obtener una mejor resolución de frecuencia y una estimación menos sesgada, use un búfer más largo, pero no tanto como para que la señal ya no sea aproximadamente estacionaria. Aquí hay un ejemplo en Python:
from pylab import *
import wave
fs = 44100.0 # sample rate
K = 3 # number of windows
L = 8192 # 1st pass window overlap, 50%
M = 16384 # 1st pass window length
N = 32768 # 1st pass DFT lenth: acyclic correlation
# load a sample of guitar playing an open string 6
# with a fundamental frequency of 82.4 Hz (in theory),
# but this sample is actually at about 81.97 Hz
g = fromstring(wave.open(''dist_gtr_6.wav'').readframes(-1),
dtype=''int16'')
g = g / float64(max(abs(g))) # normalize to +/- 1.0
mi = len(g) / 4 # start index
def welch(x, w, L, N):
# Welch''s method
M = len(w)
K = (len(x) - L) / (M - L)
Xsq = zeros(N/2+1) # len(N-point rfft) = N/2+1
for k in range(K):
m = k * ( M - L)
xt = w * x[m:m+M]
# use rfft for efficiency (assumes x is real-valued)
Xsq = Xsq + abs(rfft(xt, N)) ** 2
Xsq = Xsq / K
Wsq = abs(rfft(w, N)) ** 2
bias = irfft(Wsq) # for unbiasing Rxx and Sxx
p = dot(x,x) / len(x) # avg power, used as a check
return Xsq, bias, p
# first pass: acyclic autocorrelation
x = g[mi:mi + K*M - (K-1)*L] # len(x) = 32768
w = hamming(M) # hamming[m] = 0.54 - 0.46*cos(2*pi*m/M)
# reduces the side lobes in DFT
Xsq, bias, p = welch(x, w, L, N)
Rxx = irfft(Xsq) # acyclic autocorrelation
Rxx = Rxx / bias # unbias (bias is tapered)
mp = argmax(Rxx[28:561]) + 28 # index of 1st peak in 28 to 560
# 2nd pass: cyclic autocorrelation
N = M = L - (L % mp) # window an integer number of periods
# shortened to ~8192 for stationarity
x = g[mi:mi+K*M] # data for K windows
w = ones(M); L = 0 # rectangular, non-overlaping
Xsq, bias, p = welch(x, w, L, N)
Rxx = irfft(Xsq) # cyclic autocorrelation
Rxx = Rxx / bias # unbias (bias is constant)
mp = argmax(Rxx[28:561]) + 28 # index of 1st peak in 28 to 560
Sxx = Xsq / bias[0]
Sxx[1:-1] = 2 * Sxx[1:-1] # fold the freq axis
Sxx = Sxx / N # normalize S for avg power
n0 = N / mp
np = argmax(Sxx[n0-2:n0+3]) + n0-2 # bin of the nearest peak power
# check
print "/nAverage Power"
print " p:", p
print "Rxx:", Rxx[0] # should equal dot product, p
print "Sxx:", sum(Sxx), ''/n'' # should equal Rxx[0]
figure().subplots_adjust(hspace=0.5)
subplot2grid((2,1), (0,0))
title(''Autocorrelation, R$_{xx}$''); xlabel(''Lags'')
mr = r_[:3 * mp]
plot(Rxx[mr]); plot(mp, Rxx[mp], ''ro'')
xticks(mp/2 * r_[1:6])
grid(); axis(''tight''); ylim(1.25*min(Rxx), 1.25*max(Rxx))
subplot2grid((2,1), (1,0))
title(''Power Spectral Density, S$_{xx}$''); xlabel(''Frequency (Hz)'')
fr = r_[:5 * np]; f = fs * fr / N;
vlines(f, 0, Sxx[fr], colors=''b'', linewidth=2)
xticks((fs * np/N * r_[1:5]).round(3))
grid(); axis(''tight''); ylim(0,1.25*max(Sxx[fr]))
show()
Salida:
Average Power
p: 0.0410611012542
Rxx: 0.0410611012542
Sxx: 0.0410611012542
El retraso máximo es de 538, que es 44100/538 = 81.97 Hz. El DFT acíclico de primer paso muestra el fundamental en el contenedor 61, que es 82.10 +/- 0.67 Hz. El segundo paso utiliza una longitud de ventana de 538 * 15 = 8070, por lo que las frecuencias DFT incluyen el período fundamental y los armónicos de la cadena. Esto permite una autocorrelación cíclica parcial para una mejor estimación de PSD con una menor propagación armónica (es decir, la correlación puede envolver alrededor de la ventana periódicamente).
Edición: Actualizado para usar el método de Welch para estimar la autocorrelación. La superposición de las ventanas compensa la ventana de Hamming. También calculo el sesgo cónico de la ventana de Hamming para desvincular la autocorrelación.
Edición: Se agregó una segunda pasada con correlación cíclica para limpiar la densidad espectral de potencia. Este paso usa 3 ventanas rectangulares no superpuestas, longitud 538 * 15 = 8070 (lo suficientemente corta como para ser casi estacionaria). El sesgo para la correlación cíclica es una constante, en lugar del sesgo cónico de la ventana de Hamming.