Python中提取specgram频谱图高幅值对应时间与频率点的方法

首先明确ax.specgram()返回值的对应关系：

• Pxx：形状为(M, N)的二维数组，M为频率点总数，N为时间窗总数，每个元素对应特定频率、时间点的功率谱密度幅值
• freqs：长度为M的一维数组，存储频率轴各索引对应的实际频率值
• bins：长度为N的一维数组，存储时间轴各索引对应的实际时间点值

具体实现步骤 #

• 第一步：设定高幅值判定阈值
  你可以根据实际场景选择阈值，通用方案是取功率谱密度的高分位数，比如取幅值前5%作为高幅值区域的判定标准：

import numpy as np
# 95分位数表示大于该值的元素占所有元素的5%
threshold = np.percentile(Pxx, 95)

如果需要更精准的阈值，也可以直接观察Pxx的数值分布手动指定。

• 第二步：筛选高幅值区域的索引
  通过numpy的where方法筛选出符合阈值要求的坐标索引：

# 返回值第一个数组为频率维度的索引，第二个数组为时间维度的索引
high_amp_indexes = np.where(Pxx > threshold)

• 第三步：映射索引到实际频率、时间值
  用拿到的索引去对应取值即可得到一一匹配的高幅值频率、时间点：

# 所有高幅值点对应的频率值
high_freqs = freqs[high_amp_indexes[0]]
# 所有高幅值点对应的时间值
high_times = bins[high_amp_indexes[1]]

如果需要提取连续的高幅值区域（避免离散点过多），可以借助scipy的连通域分析工具处理：

from scipy.ndimage import label
# 生成高幅值掩码
mask = (Pxx > threshold).astype(int)
# 标记连通域，num_features为连通域总数
labeled_mask, num_features = label(mask)
# 后续可按连通域的标签单独提取每个连续高幅值区域的时间、频率范围

拿到对应的时间、频率点后，你可以按时间维度分组，将每个时间点对应的高幅值频率堆叠展示，即可实现你需要的瀑布图效果。

内容的提问来源于stack exchange，提问作者Anita