Изготовление спектрограммы из микрофона

Ответ 1

Во-первых, processBlock внимание, что ваш код отображает до 100 спектрограмм (если processBlock вызывается несколько раз) друг над другом, и вы видите только последнюю. Вы можете исправить это. Кроме того, я предполагаю, что вы знаете, почему вы хотите работать с 30 мс аудио-записями. Лично я не могу придумать практического применения, где 30 мс, записанные микрофоном ноутбука, могли бы дать интересные идеи. Это зависит от того, что вы записываете и как вы запускаете запись, но эта проблема имеет отношение к реальному вопросу.

В противном случае код работает отлично. С помощью нескольких небольших изменений в функции processBlock, применяя некоторые processBlock знания, вы можете получить информативные и эстетические спектрограммы.

Итак, давайте поговорим о реальных спектрограммах. Я возьму вывод SoX в качестве ссылки. В примечании к цветовой полосе указано, что это dBFS ¹ что является логарифмической мерой (дБ - сокращение от децибела). Итак, сначала преобразуем спектрограмму в дБ:

    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)

Большой! Теперь мы получили относительно сбалансированное разрешение по обеим осям - но подождите! Почему результат такой пиксельный? На самом деле, это вся информация, которая есть в коротком временном окне 30 мс. Есть только так много способов, которыми 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)

Я лично предпочитаю 3-е решение, но вам нужно будет найти свой собственный предпочтительный компромисс времени/частоты.

Наконец, позвольте использовать цветовую карту, которая больше похожа на SoX:

plt.pcolormesh(t, f, dBS, cmap='inferno')

Краткий комментарий к следующей строке:

THRESHOLD = 40 # dB

Порог сравнивается со среднеквадратичным значением входного сигнала, которое измеряется не в дБ, а в единицах необработанной амплитуды.

¹ Видимо, FS - это сокращение от полной шкалы. dBFS означает, что показатель дБ относится к максимальному диапазону. 0 дБ - самый громкий сигнал, возможный в текущем представлении, поэтому фактические значения должны быть <= 0 дБ.

Ответ 2

ОБНОВЛЕНИЕ, чтобы сделать мой ответ более ясным и, надеюсь, дополню отличное объяснение от @kazemakase, я нашел три вещи, которые, я надеюсь, помогут:

• Использовать LogNorm:

  plt.pcolormesh(t, f, Sxx, cmap='RdBu', norm=LogNorm(vmin=Sxx.min(), vmax=Sxx.max()))
  

• использовать метод numpy fromstring

Выключает вычисление RMS не будет работать с этим методом, поскольку данные ограничены тип данных длины, а переполнения становятся отрицательными: 507 * 507 = -5095.

3. используйте colorbar(), поскольку все становится легче, когда вы можете увидеть масштаб

   plt.colorbar()
   

Исходный ответ:

У меня получился приличный результат, который воспроизводит частоту 10 кГц в вашем коде только с несколькими изменениями:

• импортировать LogNorm

  from matplotlib.colors import LogNorm
  

• Используйте LogNorm в сетке

  plt.pcolormesh(t, f, Sxx, cmap='RdBu', norm=LogNorm(vmin=Sxx.min(), vmax=Sxx.max()))
  

Это дало мне:

Вам также может потребоваться вызвать plt.close() после savefig, и я думаю, что чтение потока требует некоторой работы, поскольку более поздние изображения отбрасывали первую четверть звука.

Id также рекомендует plt.colorbar(), чтобы вы могли увидеть масштаб, который он закончил, используя

ОБНОВЛЕНИЕ: видя, как кто-то потратил время на downvote

Вот мой код для рабочей версии спектрограммы. Он захватывает пять секунд аудио и записывает их в файл спецификации и аудиофайл, чтобы вы могли сравнивать. Theres все еще много, чтобы улучшить, и его вряд ли оптимизировали: я уверен, что его падающие куски из-за времени, чтобы написать аудио и спецификации файлов. Лучшим подходом было бы использовать неблокирующий обратный вызов, и я мог бы сделать это позже

Основное отличие от исходного кода заключалось в изменении для получения данных в правильном формате для 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()

Ответ 3

Я думаю, проблема в том, что вы пытаетесь создать спектрограмму 30-миллисекундного аудиоблока, который настолько короткий, что вы можете считать сигнал стационарным.
Спектрограмма фактически является STFT, и вы можете найти это также в документации Scipy:

scipy.signal.spectrogram (x, fs = 1.0, window = ('tukey', 0.25), nperseg = нет, noverlap = нет, nfft = нет, detrend = "константа", return_onesided = True, масштабирование = "плотность", ось = -1, mode = 'psd')

Вычислить спектрограмму с последовательными преобразованиями Фурье.

Спектрограммы можно использовать как способ визуализации изменения частотного содержания нестационарных сигналов во времени.

На первом рисунке у вас есть четыре среза, которые являются результатом четырех последовательных fft на вашем сигнальном блоке, с некоторыми окнами и перекрытиями. Второй рисунок имеет уникальный срез, но он зависит от параметров спектрограммы, которые вы использовали.
Дело в том, что вы хотите сделать с этим сигналом. Какова цель алгоритма?

Ответ 4

Я не уверен, что работа непосредственно в Python - лучший способ для обработки звука и наиболее точно с FFT... [на мой взгляд, используя cython отображаются как обязательство в обработке звука с помощью python]

Вы оцениваете возможность привязки любого внешнего метода FFT (используя fftw, например) и продолжайте использовать python только для отправки задания к внешнему методу и обновить результат изображения?

Вы можете найти некоторую информацию относительно оптимизации FFT в python здесь, а также можете посмотреть scipy Реализация FFT.