python - 无法弄清楚如何计算 .wav 文件的光谱属性

Question

我试图计算出声音剪辑的平均频率以及 Q25 和 Q75 值，但遇到了问题(主要是由于我缺乏数学和 DSP 知识)。

我要离开 this答案，并且在将该答案中的代码与读取 .wav 文件相结合时遇到问题。

这是我用来记录的代码...

def record_sample(file):
    # Audio Recording
    CHUNK = 1024
    FORMAT = pyaudio.paInt16
    CHANNELS = 2
    RATE = 44100
    RECORD_SECONDS = 5
    pa = pyaudio.PyAudio()

    # Record sample.
    stream = pa.open(format=FORMAT, channels=CHANNELS, rate=RATE, input=True, frames_per_buffer=CHUNK)
    frames = []
    for _ in range(int(RATE / CHUNK * RECORD_SECONDS)):
        data = stream.read(CHUNK)
        frames.append(data)
    stream.stop_stream()
    stream.close()

    # Save to wave file.
    wf = wave.open(file, "wb")
    wf.setnchannels(CHANNELS)
    wf.setsampwidth(pa.get_sample_size(FORMAT))
    wf.setframerate(RATE)
    wf.writeframes(b''.join(frames))
    wf.close()

一切正常。这是我用于计算平均频率、Q25 和 Q75 的代码...

def spectral_properties(file):
    # Note: scipy.io.wavfile.read
    fs, data = wavfile.read(file)
    spec = np.abs(np.fft.rfft(data))
    freq = np.fft.rfftfreq(len(data), d=1 / fs)
    spec = np.abs(spec)
    amp = spec / spec.sum()
    amp_cumsum = np.cumsum(amp)
    Q25 = freq[len(amp_cumsum[amp_cumsum <= 0.25]) + 1]
    Q75 = freq[len(amp_cumsum[amp_cumsum <= 0.75]) + 1]
    print((freq * amp).sum(), Q25, Q75)

以及它产生的错误......

File "/home/horner/workspace/school/ML/machine-learning-project-mdx97/program/audio.py", line 65, in spectral_properties
    Q75 = freq[len(amp_cumsum[amp_cumsum <= 0.75]) + 1]
IndexError: index 298981 is out of bounds for axis 0 with size 110081

Answer 1

请注意，您有 2 个 channel ，这意味着您将获得二维数据。您当前的版本只是在某些操作期间展平数组，这就是它看起来有太多元素的原因。

有两种方法可以解决这个问题。第一种是仅使用其中一个 channel :

def spectral_properties(filename):

    fs, data = wavfile.read(filename)

    # use the first channel only
    if data.ndim > 1:
        data = data[:, 0]

    spec = np.abs(np.fft.rfft(data))
    freq = np.fft.rfftfreq(len(data), d=1/fs)

    assert len(spec) == len(freq)

    amp = spec / spec.sum()
    amp_cumsum = amp.cumsum()

    assert len(amp_cumsum) == len(freq)

    q25 = freq[len(amp_cumsum[amp_cumsum < 0.25])]
    q75 = freq[len(amp_cumsum[amp_cumsum < 0.75])]

    return (freq * amp).sum(), q25, q75


avg, q25, q75 = spectral_properties('foobar.wav')
print(avg, q25, q75)

第二个是保留 channel 并告诉 numpy 函数它们应该沿哪个轴操作。这也意味着计算四分位数变得不再那么简单，因为您需要为每个 channel 单独查找它们，但由于 Python 的列表推导式，它看起来几乎和以前一样简单:

def spectral_properties(filename):

    fs, data = wavfile.read(filename)

    # determine number of channels
    num_channels = data.shape[1]

    spec = np.abs(np.fft.rfft(data, axis=0))
    freq = np.fft.rfftfreq(len(data), d=1/fs)

    assert len(spec) == len(freq)

    amp = spec / spec.sum(axis=0)
    amp_cumsum = amp.cumsum(axis=0)

    assert len(amp_cumsum) == len(freq)

    q25 = [freq[len(amp_cumsum[:,j][amp_cumsum[:,j] < 0.25])] for j in range(num_channels)]
    q75 = [freq[len(amp_cumsum[:,j][amp_cumsum[:,j] < 0.75])] for j in range(num_channels)]

    return (freq[:,np.newaxis] * amp).sum(axis=0), q25, q75


avg, q25, q75 = spectral_properties('foobar.wav')
print(avg, q25, q75)

请注意，+ 1 在查找四分位数的原始表达式中存在问题。考虑除最后一个值之外的所有值都小于 0.25。因此，对于 n - 1 元素来说，不等式成立。您添加 1，因此您得到 n。但是对于长度为 n 的 freq 数组来说，n 的索引太高了。

此外，我怀疑您可能想要对 spec 进行平方，而不是保留其大小。

更新:

您可能还想使用 searchsorted 来查找四分位数，这应该更快、更容易阅读:

q25 = freq[np.searchsorted(amp_cumsum, 0.25)]
q75 = freq[np.searchsorted(amp_cumsum, 0.75)]

还有:

q25 = [freq[np.searchsorted(amp_cumsum[:,j], 0.25)] for j in range(num_channels)]
q75 = [freq[np.searchsorted(amp_cumsum[:,j], 0.75)] for j in range(num_channels)]