单吉他音符的谐波乘积谱 Python

Harmonic product spectrum for single guitar note Python

我正在尝试检测用吉他弹奏的 B3 音符的音高。可以找到音频 here

这是频谱图:

如您所见,可以看出基本音高约为 250Hz,对应于 B3 音符。

它还包含大量谐波,这就是我选择使用 here 中的 HPS 的原因。我正在使用此代码检测音高:

def freq_from_hps(signal, fs):
    """Estimate frequency using harmonic product spectrum
    Low frequency noise piles up and overwhelms the desired peaks
    """
    N = len(signal)
    signal -= mean(signal)  # Remove DC offset

    # Compute Fourier transform of windowed signal
    windowed = signal * kaiser(N, 100)

    # Get spectrum
    X = log(abs(rfft(windowed)))

    # Downsample sum logs of spectra instead of multiplying
    hps = copy(X)
    for h in arange(2, 9): # TODO: choose a smarter upper limit
        dec = decimate(X, h)
        hps[:len(dec)] += dec

    # Find the peak and interpolate to get a more accurate peak
    i_peak = argmax(hps[:len(dec)])
    i_interp = parabolic(hps, i_peak)[0]

    # Convert to equivalent frequency
    return fs * i_interp / N  # Hz

我的采样率为 40000。但是,我得到的结果不是接近 250Hz(B3 注),而是 0.66Hz。这怎么可能?

我也尝试使用同一个 repo 中的自相关方法,但我也得到了像 10000Hz 这样的糟糕结果。

多亏了一个答案,我明白我必须应用一个滤波器来去除信号中的低频。我怎么做?是否有多种方法可以做到这一点,推荐使用哪一种?

状态更新:

答案建议的高通滤波器有效。如果我将答案中的函数应用于我的音频信号,它会正确显示大约 245Hz。但是,我想过滤整个信号,而不仅仅是其中的一部分。一个音符可能位于信号的中间,或者一个信号包含多个音符(我知道一个解决方案是开始检测,但我很想知道为什么这不起作用)。这就是为什么我将代码编辑为 return filtered_audio.

问题是如果我这样做,即使噪声已被正确移除(见屏幕截图)。结果我得到 0.05。

根据频谱图中谐波之间的距离,我估计音调约为 150-200 Hz。那么,为什么音调检测算法没有检测到我们可以在频谱图中肉眼看到的音调呢?我有几个猜测:

音符只持续几秒钟。开头有一个10个以上泛音的漂亮泛音栈!这些很快消失,甚至在 5 秒后就看不见了。如果您尝试估计整个信号的音高,您的估计可能会被 5-12 秒内声音的 "pitch" 污染。 尝试只计算前 1-2 秒的音调。

低频噪声太多。在频谱图中,您可以看到0到64 Hz之间有很多功率。这不是谐波的一部分,因此您可以尝试使用 高通滤波器.

将其移除

下面是一些完成这项工作的代码:

import numpy as np
from scipy.io import wavfile
from scipy import signal
import matplotlib.pyplot as plt

from frequency_estimator import freq_from_hps
# downloaded from https://github.com/endolith/waveform-analyzer/

filename = 'Vocaroo_s1KZzNZLtg3c.wav'
# downloaded from http://vocaroo.com/i/s1KZzNZLtg3c

# Parameters
time_start = 0  # seconds
time_end = 1  # seconds
filter_stop_freq = 70  # Hz
filter_pass_freq = 100  # Hz
filter_order = 1001

# Load data
fs, audio = wavfile.read(filename)
audio = audio.astype(float)

# High-pass filter
nyquist_rate = fs / 2.
desired = (0, 0, 1, 1)
bands = (0, filter_stop_freq, filter_pass_freq, nyquist_rate)
filter_coefs = signal.firls(filter_order, bands, desired, nyq=nyquist_rate)

# Examine our high pass filter
w, h = signal.freqz(filter_coefs)
f = w / 2 / np.pi * fs  # convert radians/sample to cycles/second
plt.plot(f, 20 * np.log10(abs(h)), 'b')
plt.ylabel('Amplitude [dB]', color='b')
plt.xlabel('Frequency [Hz]')
plt.xlim((0, 300))

# Apply high-pass filter
filtered_audio = signal.filtfilt(filter_coefs, [1], audio)

# Only analyze the audio between time_start and time_end
time_seconds = np.arange(filtered_audio.size, dtype=float) / fs
audio_to_analyze = filtered_audio[(time_seconds >= time_start) &
                                  (time_seconds <= time_end)]

fundamental_frequency = freq_from_hps(audio_to_analyze, fs)
print 'Fundamental frequency is {} Hz'.format(fundamental_frequency)