从麦克风生成频谱图
syb*_*0rg 15 python audio numpy matplotlib scipy
下面我的代码将从麦克风输入,如果音频块的平均值超过某个阈值,它将产生音频块的频谱图(长度为30 ms).以下是正常对话中生成的频谱图的样子:
从我所看到的情况来看,这看起来并不像我期望的频谱图看起来像音频及其环境.我期待更像下面的内容(转换为保留空间):
我正在录制的麦克风是我Macbook上的默认麦克风,有什么问题的建议吗？
record.py:
import pyaudio
import struct
import math
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
THRESHOLD = 40 # dB
RATE = 44100
INPUT_BLOCK_TIME = 0.03 # 30 ms
INPUT_FRAMES_PER_BLOCK = int(RATE * INPUT_BLOCK_TIME)
def get_rms(block):
return np.sqrt(np.mean(np.square(block)))
class AudioHandler(object):
def __init__(self):
self.pa = pyaudio.PyAudio()
self.stream = self.open_mic_stream()
self.threshold = THRESHOLD
self.plot_counter = 0
def stop(self):
self.stream.close()
def find_input_device(self):
device_index = None
for i in range( self.pa.get_device_count() ):
devinfo = self.pa.get_device_info_by_index(i)
print('Device %{}: %{}'.format(i, devinfo['name']))
for keyword in ['mic','input']:
if keyword in devinfo['name'].lower():
print('Found an input: device {} - {}'.format(i, devinfo['name']))
device_index = i
return device_index
if device_index == None:
print('No preferred input found; using default input device.')
return device_index
def open_mic_stream( self ):
device_index = self.find_input_device()
stream = self.pa.open( format = pyaudio.paInt16,
channels = 1,
rate = RATE,
input = True,
input_device_index = device_index,
frames_per_buffer = INPUT_FRAMES_PER_BLOCK)
return stream
def processBlock(self, snd_block):
f, t, Sxx = signal.spectrogram(snd_block, RATE)
plt.pcolormesh(t, f, Sxx)
plt.ylabel('Frequency [Hz]')
plt.xlabel('Time [sec]')
plt.savefig('data/spec{}.png'.format(self.plot_counter), bbox_inches='tight')
self.plot_counter += 1
def listen(self):
try:
raw_block = self.stream.read(INPUT_FRAMES_PER_BLOCK, exception_on_overflow = False)
count = len(raw_block) / 2
format = '%dh' % (count)
snd_block = np.array(struct.unpack(format, raw_block))
except Exception as e:
print('Error recording: {}'.format(e))
return
amplitude = get_rms(snd_block)
if amplitude > self.threshold:
self.processBlock(snd_block)
else:
pass
if __name__ == '__main__':
audio = AudioHandler()
for i in range(0,100):
audio.listen()
根据评论进行编辑:
如果我们将速率限制为16000 Hz并使用色度图的对数刻度,则这是在麦克风附近点击的输出:
这看起来对我来说仍然有些奇怪,但似乎也是朝着正确方向迈出的一步.
使用Sox并与我的程序生成的谱图进行比较:
kaz*_*ase 15
首先,观察您的代码将最多100个谱图(如果processBlock被多次调用)叠加在一起,您只能看到最后一个.你可能想解决这个问题.此外,我假设您知道为什么要使用30ms录音.就个人而言,我想不出一个实际的应用,其中笔记本电脑麦克风录制的30ms可以提供有趣的见解.它取决于您录制的内容以及触发录制的方式,但这个问题与实际问题相关.
否则代码完美无缺.只需对processBlock功能进行一些小的更改,应用一些背景知识,您就可以获得信息和美学光谱图.
那么让我们来谈谈实际的光谱图.我将把SoX输出作为参考.colorbar注释表明它是dBFS.我不知道FS 1代表什么,但dB代表Decibel,这是一个对数刻度.那么,让我们首先将频谱图转换为dB:
f, t, Sxx = signal.spectrogram(snd_block, RATE)
dBS = 10 * np.log10(Sxx) # convert to dB
plt.pcolormesh(t, f, dBS)
这改善了色标.现在我们看到之前隐藏的较高频段的噪声.接下来,让我们来解决时间问题.频谱图将信号分成段(默认长度为256)并计算每个段的频谱.这意味着我们具有出色的频率分辨率但是时间分辨率非常差,因为只有少数这样的段适合信号窗口(大约1300个样本长).在时间和频率分辨率之间总是需要权衡.这与不确定性原则有关.因此,让我们通过将信号分成更短的段来交换一些频率分辨率以获得时间分辨率:
f, t, Sxx = signal.spectrogram(snd_block, RATE, nperseg=64)
大!现在我们在两个轴上都获得了相对平衡的分辨率 - 但等等!为什么结果如此像素化？!实际上,这是短30ms时间窗口中的所有信息.只有很多方法可以在两个维度上分布1300个样本.但是,我们可以作弊并使用更高的FFT结果和重叠段.这样可以使结果更平滑,但不提供其他信息:
f, t, Sxx = signal.spectrogram(snd_block, RATE, nperseg=64, nfft=256, noverlap=60)
看到漂亮的光谱干涉图案.(这些模式取决于所使用的窗口函数,但不要在这里详细说明.请参阅window光谱图函数的参数来使用这些.)结果看起来不错,但实际上不包含比上一个图像更多的信息.
为了使结果更加SoX-lixe观察到SoX频谱图在时间轴上相当模糊.您可以通过使用原始的低时间分辨率(长段)获得此效果,但让它们重叠以获得平滑度:
f, t, Sxx = signal.spectrogram(snd_block, RATE, noverlap=250)
我个人更喜欢第三种解决方案,但您需要找到自己喜欢的时间/频率权衡.
最后,让我们使用更像SoX的色彩映射:
plt.pcolormesh(t, f, dBS, cmap='inferno')
以下简短评论:
THRESHOLD = 40 # dB
将阈值与输入信号的RMS进行比较,输入信号的RMS 不是以dB为单位测量,而是以原始幅度单位测量.
1显然FS是满量程的缩写.dBFS表示dB度量相对于最大范围.0 dB是当前表示中可能的最大信号,因此实际值必须<= 0 dB.
- 我不能不投票给可以解释时间/频率权衡的人 (2认同)
更新以使我的答案更清楚,并希望赞美@kazemakase的优秀解释,我发现了三件我希望会有所帮助的事情:
使用LogNorm:
Run Code Online (Sandbox Code Playgroud)plt.pcolormesh(t, f, Sxx, cmap='RdBu', norm=LogNorm(vmin=Sxx.min(), vmax=Sxx.max()))使用numpy的fromstring方法
事实证明RMS计算不适用于此方法,因为数据受约束长度数据类型且溢出变为负数:即507*507 = -5095.
当你可以看到比例时,使用colorbar()作为eveything变得更容易
Run Code Online (Sandbox Code Playgroud)plt.colorbar()
原答案:
我得到了一个不错的结果,在你的代码中播放10kHz的频率,只有几个改动:
导入LogNorm
Run Code Online (Sandbox Code Playgroud)from matplotlib.colors import LogNorm在网格中使用LogNorm
Run Code Online (Sandbox Code Playgroud)plt.pcolormesh(t, f, Sxx, cmap='RdBu', norm=LogNorm(vmin=Sxx.min(), vmax=Sxx.max()))
这给了我:
你可能还需要在savefig之后调用plt.close(),我认为流读取需要一些工作,因为后来的图像正在丢失声音的第一个四分之一.
我也建议plt.colorbar()你这样看你最终使用的比例
更新:看到有人花时间进行投票
下面是我对光谱图的工作版本的代码.它捕获五秒钟的音频并将它们写入spec文件和音频文件,以便您进行比较.Theres stilla很多要改进并且几乎没有优化:我确定它的丢弃块是因为有时间写音频和spec文件.更好的方法是使用非阻塞回调,稍后我可能会这样做
与原始代码的主要区别在于为numpy获取正确格式的数据的更改:
np.fromstring(raw_block,dtype=np.int16)
代替
struct.unpack(format, raw_block)
一旦我尝试使用以下方法将音频写入文件,这就成为一个主要问题:
scipy.io.wavfile.write('data/audio{}.wav'.format(self.plot_counter),RATE,snd_block)
这是一个很好的音乐,鼓很明显:
代码:
import pyaudio
import struct
import math
import numpy as np
from scipy import signal
import matplotlib.pyplot as plt
from matplotlib.colors import LogNorm
import time
from scipy.io.wavfile import write
THRESHOLD = 0 # dB
RATE = 44100
INPUT_BLOCK_TIME = 1 # 30 ms
INPUT_FRAMES_PER_BLOCK = int(RATE * INPUT_BLOCK_TIME)
INPUT_FRAMES_PER_BLOCK_BUFFER = int(RATE * INPUT_BLOCK_TIME)
def get_rms(block):
return np.sqrt(np.mean(np.square(block)))
class AudioHandler(object):
def __init__(self):
self.pa = pyaudio.PyAudio()
self.stream = self.open_mic_stream()
self.threshold = THRESHOLD
self.plot_counter = 0
def stop(self):
self.stream.close()
def find_input_device(self):
device_index = None
for i in range( self.pa.get_device_count() ):
devinfo = self.pa.get_device_info_by_index(i)
print('Device %{}: %{}'.format(i, devinfo['name']))
for keyword in ['mic','input']:
if keyword in devinfo['name'].lower():
print('Found an input: device {} - {}'.format(i, devinfo['name']))
device_index = i
return device_index
if device_index == None:
print('No preferred input found; using default input device.')
return device_index
def open_mic_stream( self ):
device_index = self.find_input_device()
stream = self.pa.open( format = self.pa.get_format_from_width(2,False),
channels = 1,
rate = RATE,
input = True,
input_device_index = device_index)
stream.start_stream()
return stream
def processBlock(self, snd_block):
f, t, Sxx = signal.spectrogram(snd_block, RATE)
zmin = Sxx.min()
zmax = Sxx.max()
plt.pcolormesh(t, f, Sxx, cmap='RdBu', norm=LogNorm(vmin=zmin, vmax=zmax))
plt.ylabel('Frequency [Hz]')
plt.xlabel('Time [sec]')
plt.axis([t.min(), t.max(), f.min(), f.max()])
plt.colorbar()
plt.savefig('data/spec{}.png'.format(self.plot_counter), bbox_inches='tight')
plt.close()
write('data/audio{}.wav'.format(self.plot_counter),RATE,snd_block)
self.plot_counter += 1
def listen(self):
try:
print "start", self.stream.is_active(), self.stream.is_stopped()
#raw_block = self.stream.read(INPUT_FRAMES_PER_BLOCK, exception_on_overflow = False)
total = 0
t_snd_block = []
while total < INPUT_FRAMES_PER_BLOCK:
while self.stream.get_read_available() <= 0:
print 'waiting'
time.sleep(0.01)
while self.stream.get_read_available() > 0 and total < INPUT_FRAMES_PER_BLOCK:
raw_block = self.stream.read(self.stream.get_read_available(), exception_on_overflow = False)
count = len(raw_block) / 2
total = total + count
print "done", total,count
format = '%dh' % (count)
t_snd_block.append(np.fromstring(raw_block,dtype=np.int16))
snd_block = np.hstack(t_snd_block)
except Exception as e:
print('Error recording: {}'.format(e))
return
self.processBlock(snd_block)
if __name__ == '__main__':
audio = AudioHandler()
for i in range(0,5):
audio.listen()
我认为问题是你正在尝试进行30ms音频块的频谱图,这个频谱很短,你可以认为信号是静止的.
谱图实际上是STFT,您也可以在Scipy文档中找到它:
scipy.signal.spectrogram(x,fs = 1.0,window =('tukey',0.25),nperseg = None,noverlap = None,nfft = None,detrend ='constant',return_onesided = True,scaling ='density', axis = -1,mode ='psd')
使用连续傅立叶变换计算频谱图.
频谱图可用作可视化非平稳信号频率内容随时间变化的方式.
在第一个图中,您有四个切片,它们是信号块上四个连续fft的结果,有一些窗口和重叠.第二个图有一个独特的切片,但它取决于您使用的谱图参数.
关键是你想用这个信号做什么.算法的目的是什么？