Python实现matlab的spectrogram函数

Python实现matlab的spectrogram函数

关于matlab的spectrogram函数的详细介绍可参考这篇博客

matlab代码：

s = spectrogram(x,window,noverlap,nfft)
%%
x---输入信号的向量。默认情况下，即没有后续输入参数，x将被分成8段分别做变换处理，如果x不能被平分成8段，则会做截断处理。默认情况下，其他参数的默认值为：window---窗函数，默认为nfft长度的海明窗Hamming；noverlap---每一段的重叠样本数，默认值是在各段之间产生50%的重叠；nfft---做FFT变换的长度，默认为256和大于每段长度的最小2次幂之间的最大值。另外，此参数除了使用一个常量外，还可以指定一个频率向量F；fs---采样频率，默认值归一化频率。
Window---窗函数，如果window为一个整数，x将被分成window段，每段使用Hamming窗函数加窗。如果window是一个向量，x将被分成length(window)段，每一段使用window向量指定的窗函数加窗。所以如果想获取specgram函数的功能，只需指定一个256长度的Hann窗。
noverlap---各段之间重叠的采样点数。它必须为一个小于window或length(window)的整数。其意思为两个相邻窗不是尾接着头的，而是两个窗有交集，有重叠的部分。
nfft---计算离散傅里叶变换的点数。它需要为标量。
S---输入信号x的短时傅里叶变换。它的每一列包含一个短期局部时间的频率成分估计，时间沿列增加，频率沿行增加。如果x是长度为Nx的复信号，则S为nfft行k列的复矩阵，其中k取决于window，如果window为一个标量，则k = fix((Nx-noverlap)/(window-noverlap))；如果window为向量，则k = fix((Nx-noverlap)/(length(window)-noverlap))。对于实信号x，如果nfft为偶数，则S的行数为(nfft/2+1)，如果nfft为奇数，则行数为(nfft+1)/2，列数同上。
%%

对应的python代码：

from scipy.fftpack import fft
from scipy import signal
import numpy as np

def Spectrogram(x, window, noverlap, nfft):
    
    len_x = len(x)
    step = window - noverlap
    nn = nfft//2+1
    num_win = int(np.floor((len_x - noverlap)/ (window - noverlap)))
    spectrogram_data = []
    # Hamming 窗，matlab函数spectogram默认加该窗
    win = signal.hamming(window)  
    for i in range(num_win):
        subdata = x[i*step : i*step + window]
        F = fft(subdata * win, n=nfft)
        spectrogram_data.append(F[:nn])
    spectogram_data = np.array(spectogram_data)
    # 为了使得输出矩阵行列与matlab结果一致，进行转置
    spectogram_data = spectogram_data.T
    return spectogram_data