Вычислить БПФ с определенной длиной БПФ, размером блока и окнами из волнового файла

Question

У меня есть 32-битный файл .wav с плавающей запятой, записанный в течение 5 секунд с частотой дискретизации 48 кГц. И я хочу получить БПФ для полных 5 секунд с примененным окном Хеннинга, без перекрытия и с длиной БПФ 8192. При получении БПФ я понимаю, что делаю ошибку в вычислении БПФ.

Для этого я использовал следующий код:

      # Imports
import matplotlib.pyplot as plt
import numpy as np
import math as math
import soundfile as sf  
import scipy.fftpack as fftpack
import scipy.io.wavfile as wf


def dbfft(x, fs, win=None, ref=1):
    """
    Calculate spectrum in dB scale
    Args:
        x: input signal
        fs: sampling frequency
        win: vector containing window samples (same length as x).
             If not provided, then rectangular window is used by default.
        ref: reference value used for dBFS . 32768 for int16 and 1 for float

    Returns:
        freq: frequency vector
        s_db: spectrum in dB scale
    """
    N = len(x)  # Length of input sequence
    if win is None:
        win = np.ones(1, N)
    if len(x) != len(win):
            raise ValueError('Signal and window must be of the same length')
    x = x * win
    # Calculate real FFT and frequency vector
    sp = np.fft.rfft(x)
    freq = np.arange((N / 2) + 1) / (float(N) / fs)
    # Scale the magnitude of FFT by window and factor of 2,because we are using half of FFT spectrum.
    s_mag = np.abs(sp) * 2 / np.sum(win)
    # Convert to dBFS using 20*log10(val/max)
    s_dbfs = 20 * np.log10(s_mag/ref)      
    return freq, s_dbfs

def main():
    # Read from wav file
    data, samplerate = sf.read('Signal Analyzer_5s_55_1_L1F4_APxRecordedAudio.wav')
    # Scaling factor
    K = 120
    # Calculation
    N = 8192
    win = np.hanning(N)
    # Frequency domain
    freq, s_dbfs = dbfft(data[0:N],samplerate, win)
    # Scale from dbFS to dB
    s_db = s_dbfs + K
    # Time domain
    Time = np.linspace(0, len(data) / samplerate, num=len(data))
    # Amp_time = 20 * np.log10 (abs(data) / 0.00002)    # ref = 1 
    Amp_time = (20*np.log10((data/50))) + K + 20*np.log10(math.sqrt(2))      # reference of sound i.e 20 log10(P rms_value/P ref)+120 dB TODO
  

    # Plot
    #fig, axes = plt.subplots(nrows=2, ncols=1)
    plt.subplot(2,1,1)
    plt.plot(Time, Amp_time)
    plt.grid('on')
    plt.minorticks_on
    plt.xlabel('Time [s]')
    plt.ylabel('Instantaneous Level [dBFS]')
    plt.xlim([0,5])
    plt.subplot(2,1,2)
    plt.plot(freq, s_db)
    plt.grid('on')
    plt.minorticks_on
    plt.xlabel('Frequency [Hz]')
    plt.ylabel('Amplitude [dB]')
    plt.xscale('log')
    plt.xlim([10,10000])
    plt.show()
    
if __name__ == "__main__":
    main()

Во время кода я вижу, что делаю БПФ только для первых 8192 отсчетов и усредняю ​​БПФ для полных 240000(5 секунд) отсчетов при размере блока 8192 с окном Хэннинга. Должен ли я выполнять несколько БПФ каждые 8192 секунды в течение 5 секунд (88 БПФ) и усреднять амплитуды, чтобы получить результирующее БПФ? Есть ли эффективный способ сделать это?

Answer 1

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

Возможно, вы можете начать работу с librosa . Например, построение спектрограммы с помощью specshow .