我从WAV文件获取峰值频率。但对于录制的双通道WAV文件，它无法正常工作。

Question

我将从wav文件中获取峰值频率。

我用以下代码从wav文件获取峰值频率：

import wave
import struct
import numpy as np
import wave
import contextlib

if __name__ == '__main__':
    fname = "test.wav"
    frate = 0
    data_size = 0
    with contextlib.closing(wave.open(fname,'r')) as f:
        frate = f.getframerate()
        data_size = f.getnframes()
    wav_file = wave.open(fname, 'r')
    data = wav_file.readframes(data_size)
    data_size = data_size * wav_file.getnchannels()
    print wav_file.getparams()
    wav_file.close()
    data = struct.unpack('{n}h'.format(n=data_size), data)
    data = np.array(data)

    w = np.fft.fft(data)
    freqs = np.fft.fftfreq(len(w))
    print(freqs.min(), freqs.max())

    # Find the peak in the coefficients
    idx = np.argmax(np.abs(w))
    freq = freqs[idx]
    freq_in_hertz = abs(freq * frate)
    print(freq_in_hertz)

我录制了一个采样率为48000、位宽为16比特、2个通道的wav文件。

在该文件中，我有一个1000Hz的正弦波。

但是脚本输出的只有500Hz。

我不知道哪里出错了。

但是对于单声道和生成的采样率为48000、位宽为16比特、2个通道的wav文件，它可以正常工作。

我使用以下脚本生成了wav文件：

import math
import wave
import struct

if __name__ == '__main__':
    # https://dev59.com/7HA65IYBdhLWcg3wyRuk
    # http://www.sonicspot.com/guide/wavefiles.html
    freq = 1000
    data_size = 454656 * 2
    fname = "test.wav"
    frate = 48000.0
    amp = 64000.0
    nchannels = 2
    sampwidth = 2
    framerate = int(frate)
    nframes = data_size
    comptype = "NONE"
    compname = "not compressed"
    data = [math.sin(2 * math.pi * freq * (x / frate))
            for x in range(data_size)]
    wav_file = wave.open(fname, 'w')
    wav_file.setparams(
        (nchannels, sampwidth, framerate, nframes, comptype,     compname))
    for v in data:
        wav_file.writeframes(struct.pack('h', int(v * amp / 2)))
    wav_file.close()

我不知道我哪里错了。 我上传了我的wav文件到脚本生成的wav script_gen.wav，采样率为48000，2个声道，16位。 录制的wavs：2声道wav 采样率为48000，2个声道，16位1声道wav(不允许在此处发布链接，因此将在评论中发布) 采样率为48000，1个声道，16位。

我在audacity中检查了所有这些峰值频率，只显示1000Khz。

但是当我尝试使用我的脚本时，我得到了1声道wav的正确输出，而2声道wav则失败了。

更新：对于2个声道，我得到了一半的峰值频率作为输出。

我感觉自己错过了什么。 有人可以帮帮我吗？

- user6464594

---

1 channel wav - user6464594

Answer 1

为什么这么复杂？考虑以下内容。

#!/usr/bin/env python3
import numpy as np
from numpy import fft
import scipy.io.wavfile as wf
import matplotlib.pyplot as plt

sr = 44100    # sample rate
len_sig = 2   # length of resulting signal in seconds

f = 1000      # frequency in Hz

# set you time axis
t = np.linspace(0, len_sig, sr*len_sig)

# set your signal
mono_data = np.sin(2*np.pi*t*f)

# write single channel .wav file
wf.write('mono.wav', sr, mono_data)

# write two-channel .wav file 
stereo_data = np.vstack((mono_data, mono_data)).T
wf.write('stereo.wav', sr, stereo_data)

现在通过加载和分析数据来测试它。

# Load data
mono_sr, mono_data = wf.read('mono.wav')
stereo_sr, stereo_data = wf.read('stereo.wav')

# analyze the data
X_mono = fft.fft(mono_data) / len(mono_data)    # remember to normalize your amplitudes

# Remember that half of energy of the signal is distributed over the 
# positive frequencies and the other half over the negative frequencies.
# 
# Commonly you want see a magnitude spectrum. That means, we ignore the phases. Hence, we
# simply multiply the spectrum by 2 and consider ONLY the first half of it.
freq_nq = len(X_mono) // 2
X_mono = abs(X_mono[:freq_nq]) * 2
freqs_mono = fft.fftfreq(len(mono_data), 1/mono_sr)[:freq_nq]

# in order the analyze a stereo signal you first have to add both channels
sum_stereo = stereo_data.sum(axis=1) / 2

# and now the same way as above
freq_nq = len(sum_stereo) // 2
X_stereo= abs(fft.fft(sum_stereo))[:freq_nq] / len(stereo_data) * 2
freqs_stereo = fft.fftfreq(len(stereo_data), 1/stereo_sr)[:freq_nq]

峰值拾取：

freqs_mono[np.argmax(X_mono)]        # == 1000.0
freqs_stereo[np.argmax(X_stereo)]    # == 1000.0

绘制结果：

fig, (ax1, ax2) = plt.subplots(2, figsize=(10,5), sharex=True, sharey=True)
ax1.set_title('mono signal')
ax1.set_xlim([0, 2000])
ax1.plot(freqs_mono, X_mono, 'b', lw=2)

ax2.set_title('stereo signal')
ax2.plot(freqs_stereo, X_stereo, 'g', lw=2)
ax2.set_xlim([0, 2000])
plt.tight_layout()
plt.show()

- MaxPowers

---

这里的“wf”和“fft”是什么意思？你导入了哪些库？ - optimus prime

@vvn 我更新了代码，包括缺失的导入语句。 - MaxPowers

@SoCelectron 我添加了最简单的选择峰值的方法。 - MaxPowers

这个错误可能是因为你正在使用Python 2.7，其中/运算符对于整数操作数执行整数除法。你可以通过在行freqs_mono = fft.fftfreq(len(mono_data), 1./mono_sr)[:freq_nq]中的1后面添加一个点来修复它。 - MaxPowers

显然，你的“data”是一维的。因此，你无法进行求和操作。 - MaxPowers

Answer 2

我认为这会帮助你走得更远。 我添加了一些额外的东西，以便与你所寻找的方式一起使用。 使用了MaxPowers的逻辑。 你需要将24位数据转换为32位，然后这也可以用于24位。

import sys
import wave
import struct
import numpy as np
import wave
import argparse

def parse_arguments():
    """Parses command line arguments."""
    parser = argparse.ArgumentParser(description='Tool to get peak frequency')
    parser.add_argument('fname', metavar='test.wav', type=str,
                        help='Path to a wav file')
    args = parser.parse_args()
    return args


def main():
    args = parse_arguments()
    fname = args.fname
    wav_file = wave.open(fname, 'r')
    frate = wav_file.getframerate()
    data_size = wav_file.getnframes()
    data = wav_file.readframes(data_size)
    nChannels = wav_file.getnchannels()
    nSample = wav_file.getsampwidth()
    data_size = data_size * nChannels * nSample
    wav_file.close()
    if nSample == 2:
        fmt = "<i2"
    else :
        fmt = "<i4"
    data = np.frombuffer(data,dtype=fmt)
    if nChannels == 2 :
        data = data.reshape(-1,nChannels)
        data = data.sum(axis=1) / 2
    # and now the same way as above as said by maxpowers
    freq_nq = len(data) // 2
    X= abs(np.fft.fft(data))[:freq_nq] / len(data) * 2
    freqs = np.fft.fftfreq(len(data), 1./frate)[:freq_nq]
    print freqs[np.argmax(X)] 

if __name__ == '__main__':
    try:
        main()
    except (Exception) as e:
        print str(e)
        sys.exit(255)

- optimus prime