本教程为脑机学习者Rose原创(转载请联系作者授权)发表于公众号：脑机接口社区(微信号：Brain_Computer).QQ交流群：903290195
倒频谱定义
倒频谱可以分析复杂频谱图上的周期结构，分离和提取在密集调频信号中的周期成分，对于具有同族谐频、异族谐频和多成分边频等复杂信号的分析非常有效。倒频谱变换是频域信号的傅立叶积分变换的再变换。时域信号经过傅立叶积分变换可转换为频率函数或功率谱密度函数，如果频谱图上呈现出复杂的周期结构而难以分辨时，对功率谱密度取对数再进行一次傅立叶积分变换，可以使周期结构呈便于识别的谱线形式。第二次傅立叶变换的平方就是倒功率谱，即“对数功率谱的功率谱”。倒功率谱的开方即称幅值倒频谱，简称倒频谱。
简言之，倒频谱分析技术是将时域振动信号的功率谱对数化，然后进行逆傅里叶变化后得到的。倒频谱的水平轴为“倒频率”的伪时间，垂直轴为对应倒频率的幅值，其计算公式为：
倒频谱python案例
实现如下：
from scipy.fftpack import fft, fftshift, ifft
from scipy.fftpack import fftfreq
import numpy as np
import matplotlib.pyplot as plt
import warnings
warnings.filterwarnings("ignore")
fs = 1000
#采样点数
num_fft = 1024
"""
生成原始信号序列
在原始信号中加上噪声
np.random.randn(t.size)
其中y1是主频为5/10/20Hz的低频信号+噪声信号;
y2是主频为50、100/200Hz的高频信号+噪声信号;
y是y1和y2的调制结果
"""
t = np.arange(0, 5, 1/fs)
y1 = 10*np.cos(2*np.pi*5*t) + 7*np.cos(2*np.pi*10*t) + 5*np.cos(2*np.pi*20*t) + np.random.randn(t.size)
y2 = 20*np.cos(2*np.pi*50*t) + 15*np.cos(2*np.pi*100*t) + 25*np.cos(2*np.pi*200*t) + np.random.randn(t.size)
y = y1*y2
plt.figure(figsize=(20, 12))
ax=plt.subplot(331)
ax.set_title('y1')
plt.plot(y1)
ax=plt.subplot(332)
ax.set_title('y2')
plt.plot(y2)
ax=plt.subplot(333)
ax.set_title('y=y1*y2')
plt.plot(y)
"""
对低频信号y1进行 FFT(Fast Fourier Transformation)快速傅里叶变换
"""
Y1 = fft(y1, num_fft)
Y1 = np.abs(Y1)
ax=plt.subplot(334)
ax.set_title('y1 fft')
plt.plot(Y1[:num_fft//2])
"""
对高频信号y2进行 FFT
"""
Y2 = fft(y2, num_fft)
Y2 = np.abs(Y2)
ax=plt.subplot(335)
ax.set_title('y2 fft')
plt.plot(Y2[:num_fft//2])
"""
对信号y进行 FFT
"""
Y = fft(y, num_fft)
Y = np.abs(Y)
ax=plt.subplot(336)
ax.set_title('y fft')
plt.plot(Y[:num_fft//2])
plt.tight_layout()
plt.show()
"""
倒频谱的定义表述为：信号→功率谱→对数→傅里叶逆变换
"""
spectrum = np.fft.fft(y, n=num_fft)
ceps = np.fft.ifft(np.log(np.abs(spectrum))).real
plt.figure(figsize=(10, 5))
plt.plot(np.abs(ceps)[:num_fft//2])
plt.title('y->spectrum->log->ifft')
plt.ylim([0, 0.2])
plt.show()
代码来源于网络，本文对代码进行注释并整理
脑机学习者Rose笔记分享,QQ交流群：903290195
更多分享，请关注公众号