从幅度或 FFT 到 dB
From Amplitude or FFT to dB
我有一个 Python 代码,它对 wav 文件执行 FFT 并绘制振幅与时间/振幅与频率图。我想从这些图中计算分贝(它们是长数组)。我不想计算精确的 dBA,我只想在计算后看到线性关系。我有分贝计,我会比较它。这是我的代码:
#!/usr/bin/env python
# -*- coding: utf-8 -*-
from __future__ import print_function
import scipy.io.wavfile as wavfile
import scipy
import scipy.fftpack
import numpy as np
from matplotlib import pyplot as plt
fs_rate, signal = wavfile.read("output.wav")
print ("Frequency sampling", fs_rate)
l_audio = len(signal.shape)
print ("Channels", l_audio)
if l_audio == 2:
signal = signal.sum(axis=1) / 2
N = signal.shape[0]
print ("Complete Samplings N", N)
secs = N / float(fs_rate)
print ("secs", secs)
Ts = 1.0/fs_rate # sampling interval in time
print ("Timestep between samples Ts", Ts)
t = scipy.arange(0, secs, Ts) # time vector as scipy arange field / numpy.ndarray
FFT = abs(scipy.fft(signal))
FFT_side = FFT[range(N//4)] # one side FFT range
freqs = scipy.fftpack.fftfreq(signal.size, t[1]-t[0])
fft_freqs = np.array(freqs)
freqs_side = freqs[range(N//4)] # one side frequency range
fft_freqs_side = np.array(freqs_side)
makespositive = signal[44100:]*(-1)
logal = np.log10(makespositive)
sn1 = np.mean(logal[1:44100])
sn2 = np.mean(logal[44100:88200])
sn3 = np.mean(logal[88200:132300])
sn4 = np.mean(logal[132300:176400])
print(sn1)
print(sn2)
print(sn3)
print(sn4)
abs(FFT_side)
for a in range(500):
FFT_side[a] = 0
plt.subplot(311)
p1 = plt.plot(t[44100:], signal[44100:], "g") # plotting the signal
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.subplot(312)
p1 = plt.plot(t[44100:], logal, "r") # plotting the signal
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.subplot(313)
p3 = plt.plot(freqs_side, abs(FFT_side), "b") # plotting the positive fft spectrum
plt.xlabel('Frequency (Hz)')
plt.ylabel('Count single-sided')
plt.show()
第一个图是幅度与时间的关系,第二个是前一个图的对数,最后一个是 FFT。
在 sn1,sn2 部分,我试图从信号中计算 dB。首先我记录日志,然后计算每秒的平均值。它没有给我一个明确的关系。这个我也试过了,没用。
import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wf
fs, signal = wf.read('output.wav') # Load the file
ref = 32768 # 0 dBFS is 32678 with an int16 signal
N = 8192
win = np.hamming(N)
x = signal[0:N] * win # Take a slice and multiply by a window
sp = np.fft.rfft(x) # Calculate real FFT
s_mag = np.abs(sp) * 2 / np.sum(win) # Scale the magnitude of FFT by window and factor of 2,
# because we are using half of FFT spectrum
s_dbfs = 20 * np.log10(s_mag / ref) # Convert to dBFS
freq = np.arange((N / 2) + 1) / (float(N) / fs) # Frequency axis
plt.plot(freq, s_dbfs)
plt.grid(True)
那么我应该执行哪些步骤? (Sum/mean 所有频率振幅然后取对数或反转,或对信号等执行)
import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wf
fs, signal = wf.read('db1.wav')
signal2 = signal[44100:]
chunk_size = 44100
num_chunk = len(signal2) // chunk_size
sn = []
for chunk in range(0, num_chunk):
sn.append(np.mean(signal2[chunk*chunk_size:(chunk+1)*chunk_size].astype(float)**2))
print(sn)
logsn = 20*np.log10(sn)
print(logsn)
输出:
[4.6057844427695475e+17, 5.0025315250895744e+17, 5.028593412665193e+17, 4.910948397471887e+17]
[353.26607217 353.98379668 354.02893044 353.82330741]
分贝计测量信号的平均功率。因此,根据您的时间信号记录,您可以使用以下公式计算平均信号功率:
chunk_size = 44100
num_chunk = len(signal) // chunk_size
sn = []
for chunk in range(0, num_chunk):
sn.append(np.mean(signal[chunk*chunk_size:(chunk+1)*chunk_size]**2))
然后以分贝为单位的相应平均信号功率简单地给出:
logsn = 10*np.log10(sn)
使用 Parseval's theorem 也可以获得频域信号的等效关系,但在您的情况下需要不必要的 FFT 计算(当您已经必须计算 FFT 时,此关系最有用用于其他目的)。
但是请注意,根据您比较的内容,可能会有一些(希望很小)差异。例如使用非线性放大器和扬声器会影响关系。同样,环境噪声会增加分贝计的测量功率。
我有一个 Python 代码,它对 wav 文件执行 FFT 并绘制振幅与时间/振幅与频率图。我想从这些图中计算分贝(它们是长数组)。我不想计算精确的 dBA,我只想在计算后看到线性关系。我有分贝计,我会比较它。这是我的代码:
#!/usr/bin/env python
# -*- coding: utf-8 -*-
from __future__ import print_function
import scipy.io.wavfile as wavfile
import scipy
import scipy.fftpack
import numpy as np
from matplotlib import pyplot as plt
fs_rate, signal = wavfile.read("output.wav")
print ("Frequency sampling", fs_rate)
l_audio = len(signal.shape)
print ("Channels", l_audio)
if l_audio == 2:
signal = signal.sum(axis=1) / 2
N = signal.shape[0]
print ("Complete Samplings N", N)
secs = N / float(fs_rate)
print ("secs", secs)
Ts = 1.0/fs_rate # sampling interval in time
print ("Timestep between samples Ts", Ts)
t = scipy.arange(0, secs, Ts) # time vector as scipy arange field / numpy.ndarray
FFT = abs(scipy.fft(signal))
FFT_side = FFT[range(N//4)] # one side FFT range
freqs = scipy.fftpack.fftfreq(signal.size, t[1]-t[0])
fft_freqs = np.array(freqs)
freqs_side = freqs[range(N//4)] # one side frequency range
fft_freqs_side = np.array(freqs_side)
makespositive = signal[44100:]*(-1)
logal = np.log10(makespositive)
sn1 = np.mean(logal[1:44100])
sn2 = np.mean(logal[44100:88200])
sn3 = np.mean(logal[88200:132300])
sn4 = np.mean(logal[132300:176400])
print(sn1)
print(sn2)
print(sn3)
print(sn4)
abs(FFT_side)
for a in range(500):
FFT_side[a] = 0
plt.subplot(311)
p1 = plt.plot(t[44100:], signal[44100:], "g") # plotting the signal
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.subplot(312)
p1 = plt.plot(t[44100:], logal, "r") # plotting the signal
plt.xlabel('Time')
plt.ylabel('Amplitude')
plt.subplot(313)
p3 = plt.plot(freqs_side, abs(FFT_side), "b") # plotting the positive fft spectrum
plt.xlabel('Frequency (Hz)')
plt.ylabel('Count single-sided')
plt.show()
第一个图是幅度与时间的关系,第二个是前一个图的对数,最后一个是 FFT。 在 sn1,sn2 部分,我试图从信号中计算 dB。首先我记录日志,然后计算每秒的平均值。它没有给我一个明确的关系。这个我也试过了,没用。
import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wf
fs, signal = wf.read('output.wav') # Load the file
ref = 32768 # 0 dBFS is 32678 with an int16 signal
N = 8192
win = np.hamming(N)
x = signal[0:N] * win # Take a slice and multiply by a window
sp = np.fft.rfft(x) # Calculate real FFT
s_mag = np.abs(sp) * 2 / np.sum(win) # Scale the magnitude of FFT by window and factor of 2,
# because we are using half of FFT spectrum
s_dbfs = 20 * np.log10(s_mag / ref) # Convert to dBFS
freq = np.arange((N / 2) + 1) / (float(N) / fs) # Frequency axis
plt.plot(freq, s_dbfs)
plt.grid(True)
那么我应该执行哪些步骤? (Sum/mean 所有频率振幅然后取对数或反转,或对信号等执行)
import numpy as np
import matplotlib.pyplot as plt
import scipy.io.wavfile as wf
fs, signal = wf.read('db1.wav')
signal2 = signal[44100:]
chunk_size = 44100
num_chunk = len(signal2) // chunk_size
sn = []
for chunk in range(0, num_chunk):
sn.append(np.mean(signal2[chunk*chunk_size:(chunk+1)*chunk_size].astype(float)**2))
print(sn)
logsn = 20*np.log10(sn)
print(logsn)
输出:
[4.6057844427695475e+17, 5.0025315250895744e+17, 5.028593412665193e+17, 4.910948397471887e+17]
[353.26607217 353.98379668 354.02893044 353.82330741]
分贝计测量信号的平均功率。因此,根据您的时间信号记录,您可以使用以下公式计算平均信号功率:
chunk_size = 44100
num_chunk = len(signal) // chunk_size
sn = []
for chunk in range(0, num_chunk):
sn.append(np.mean(signal[chunk*chunk_size:(chunk+1)*chunk_size]**2))
然后以分贝为单位的相应平均信号功率简单地给出:
logsn = 10*np.log10(sn)
使用 Parseval's theorem 也可以获得频域信号的等效关系,但在您的情况下需要不必要的 FFT 计算(当您已经必须计算 FFT 时,此关系最有用用于其他目的)。
但是请注意,根据您比较的内容,可能会有一些(希望很小)差异。例如使用非线性放大器和扬声器会影响关系。同样,环境噪声会增加分贝计的测量功率。