标签：频谱 200 Fs 浅谈 Python fft fftlw Matlab sin
浅谈快速傅里叶变化的Matlab与Python实现与区别
信号处理免不了要求频率、画频谱图，但Matlab的 fft() 函数与Python的 numpy.fft.fft() 与 scipy.fftpack.fft() 函数得到的是fft变化后的双边复数值，离画频谱图还有几句代码的距离。
基本原理不介绍了，下面直接懒人投喂，给出Matlab与Python的两个函数，直接调用即可画频谱图。
注：两种语言的fft算法是有区别的，最后细聊！
Matlab的fftlw函数
输入是信号序列、对应的时间序列、以及是否作图，输出可以得到单边归一化之后的频率与对应的振幅，通过输出可以直接画出常用的频谱图！
function [ F,M ] = fftlw( x,y,draw )
%FFTLW 快速傅里叶变化2021.10.26
%输入 x--时间 y--信号 draw--1为画频谱图，0为不画
%输出 F--频率 M--幅值
N=length(y); %采样点数
if(mod(N,2)>0)
N=N-1;
end
Fs=(N-1)/(x(N)-x(1)); %采样频率
F=(N/2:N-1)*Fs/N-Fs/2 ; %频率
y2=abs(fftshift(fft(y(1:N)))); %快速傅里叶变化
M=2*y2(N/2+1:N)/N; %归一化
M(1)=M(1)/2; %常量除以2
if draw==1 %可视化
figure
plot(F,M)
xlabel('f/HZ')
ylabel('amplitude')
title('频谱图')
end
end
Python的fftlw函数
输入与matlab的略有点不同，分别是采样频率、信号序列、是否作图，输出与matlab的函数一致。
import numpy as np
from scipy.fftpack import fft,ifft
import matplotlib.pyplot as plt
def fftlw(Fs,y,draw):
'''
Parameters
----------
Fs : 采样频率
y : 信号序列
draw :1为画频谱图，0为不画
Returns
-------
f : 频率
M : 幅值
'''
L=len(y) #采样点数
f = np.arange(int(L / 2)) * Fs / L #频率
#M = np.abs(np.fft.fft(y))*2/L #采用numpy.fft.fft()函数并归一化
M = np.abs((fft(y))) *2/L #采用scipy.fftpack.fft()函数并归一化
M = M[0:int(L / 2)] #取单边谱
M[0]=M[0]/2 #常量除以2
if draw==1: #可视化
plt.figure()
plt.rcParams['font.sans-serif']=['SimHei']
plt.rcParams['axes.unicode_minus'] = False
plt.plot(f,M)
plt.xlabel('f/HZ')
plt.ylabel('amplitude')
plt.title('频谱图')
return f,M
构造简单的信号对比两种语言fftlw效果
举个例子，构造如下信号验证所写函数的正确性：
y = 3 + t ⋅ sin ( 2 π t ⋅ 100 ) + 3 ⋅ sin ( 2 π t ⋅ 200 ) y=3+t\cdot \sin (2\pi t\cdot 100)+3\cdot \sin (2\pi t\cdot 200) y=3+t⋅sin(2πt⋅100)+3⋅sin(2πt⋅200)
其中，包括常数项，周期项和趋势项，理论上该信号的频率应该为0Hz、100Hz、200Hz（具体怎么算的补一补书知识）。在这里，我设置采样频率 fs=10000，产生10000个数据点，时域图如下：
Matlab调用fftlw函数的主函数
fs=10000; %采样频率
x=0:1/fs:(10000-1)/fs; %时间序列
y=sin(2*pi*x*100).*x+3*sin(2*pi*x*200)+3; %信号序列
figure %画时域图
plot(x,y)
title('时域图')
xlabel('t/s')
ylabel('amplitude')
[f,m]=fftlw(x,y,1); %快速傅里叶变化并画频谱图
得到的频谱图如下：发现0Hz、100Hz、200Hz处的幅值分别为3，0.5，3。0Hz与200Hz处的幅值完美对应时域中常数项与
s
i
n
(
2
π
t
⋅
200
)
sin (2\pi t\cdot 200)
sin(2πt⋅200)的系数；而100Hz项
s
i
n
(
2
π
t
⋅
200
)
sin (2\pi t\cdot 200)
sin(2πt⋅200)的系数不是常数而是
t
t
t，且时间是0-1s，该项傅里叶变化得到的是该段时间内的平均幅值，也就是0.5。
Python调用fftlw函数的主函数
直接加在def fftlw()的后文调用他就行。
Fs=10000 #采用频率
L=10000 #采样点数
t=np.arange(0,L/Fs,1/Fs) #时间序列
y=np.sin(2*np.pi*t*100)*t+3*np.sin(2*np.pi*t*200)+3 #信号序列
f,M=fftlw(Fs,y,1) #快速傅里叶变化并画频谱图
图和matlab的一模一样！但是！
采用实际的振动信号对比两种语言fftlw效果
数据来源于西储大学转子实验台振动信号，采样频率为12000Hz，现取正常状态下、转速1796 rpm轴承振动信号1000个点如下。粗略的观察，有一个低频信号大概周期为0.035 s，频率大概 29Hz。
Matlab画频谱图
Python画频谱图
python的频谱图的幅值与原始数据量级差别较大，与matlab的频谱图也毫不相关，可能是底层傅里叶变换的算法不同所致，具体哪个正确还带进一步考证！！！
标签：频谱,200,Fs,浅谈,Python,fft,fftlw,Matlab,sin 来源： https://blog.csdn.net/Lwwwwwwwl/article/details/120988992
本站声明： 1. iCode9 技术分享网（下文简称本站）提供的所有内容，仅供技术学习、探讨和分享； 2. 关于本站的所有留言、评论、转载及引用，纯属内容发起人的个人观点，与本站观点和立场无关； 3. 关于本站的所有言论和文字，纯属内容发起人的个人观点，与本站观点和立场无关； 4. 本站文章均是网友提供，不完全保证技术分享内容的完整性、准确性、时效性、风险性和版权归属；如您发现该文章侵犯了您的权益，可联系我们第一时间进行删除； 5. 本站为非盈利性的个人网站，所有内容不会用来进行牟利，也不会利用任何形式的广告来间接获益，纯粹是为了广大技术爱好者提供技术内容和技术思想的分享性交流网站。