如何获得声音文件特定频率的功率?
Posted
技术标签:
【中文标题】如何获得声音文件特定频率的功率?【英文标题】:How do I get the power at a particular frequency of a sound file? 【发布时间】:2021-02-18 08:40:39 【问题描述】:我的学位论文即将结束,我必须在特定频率 (2000Hz) 下测量水下录音(wav 文件)的声压级。所以我想出了这个代码:
''' def get_value(filename, f0, NFFT=8192, plot = False):
#Load audio
data, sampling_frequency = soundfile.read(filename)
# remove stereo
if len(data.shape)> 1:
data = data[:, 0]
# remove extra length
if len(data)>sampling_frequency:
data = data[0:sampling_frequency]
# remove DC
data = data - data.mean()
# power without filtering
total_power = 10*np.log10(np.mean(data**2))
# fft
NFFT = 4096 # number of samples in the FFT
window = np.array(1) #np.hamming(len(data))
fftdata = np.fft.fft(data / NFFT, n = NFFT)
SPL = 20 * np.log10(np.abs(fftdata)) # Sound Pressure Level [dB]
freq = np.linspace(0, sampling_frequency, NFFT) # frequency axis [Hz]
# take value at desired frequency
power_at_frequency = SPL[np.argmin(np.abs(freq-f0))]
print(power_at_frequency)
''' 但是,我大胆地检查了这个值,结果完全不同。
先谢谢了。
【问题讨论】:
【参考方案1】:如果您只对一个频率感兴趣,则无需计算 FFT,您可以简单地使用
totalEnergy = np.sum((data - np.mean(data)) ** 2)
freqEnergy = np.abs(np.sum(data * np.exp(2j * np.pi * np.arange(len(data)) * target_freq / sampling_freq)))
如果您使用 FFT 并且窗口大小不是波周期的倍数,则频率将泄漏到其他频率。为避免这种情况,您的
import numpy as np;
import matplotlib.pyplot as plt
sampling_frequency = 48000;
target_frequency = 2000.0;
ns = 1000000;
data = np.sin(2*np.pi * np.arange(ns) * target_frequency / sampling_frequency);
# power
print('a sine wave have power 0.5 ~', np.mean(data**2), 'that will be split in two ')
## Properly scaled frequency
plt.figure(figsize=(12, 5))
plt.subplot(121);
z = np.abs(np.fft.fft(data[:8192])**2) / 8192**2
print('tuned with 8192 samples', max(z), ' some power leaked in other frequencies')
plt.semilogy(np.fft.fftfreq(len(z)) * sampling_frequency, z)
plt.ylabel('power')
plt.title('some power leaked')
plt.subplot(122);
# 6000 samples = 1/8 second is multiple of 1/2000 second
z = np.abs(np.fft.fft(data[:6000])**2) / 6000**2
print('tuned with 6000 samples', max(z))
plt.semilogy(np.fft.fftfreq(len(z)) * sampling_frequency, z)
plt.xlabel('frequency')
plt.title('all power in exact two symmetric bins')
## FFT of size not multiple of 2000
print(np.sum(np.abs(np.fft.fft(data[:8192]))**2) / 8192)
【讨论】:
好的,我尝试了你的非 fft 方式,但功率值没有意义。我将功率转换为分贝单位,它们都变成了正值,而它们都应该在 -40dB 左右。 -40dB 参考哪个?你能分享你的wave文件吗?以上是关于如何获得声音文件特定频率的功率?的主要内容,如果未能解决你的问题,请参考以下文章
如何在 Linux (Ubuntu) 中输出频率 1kHz 和功率/音量 = 60%?