《Matlab如何进行离散傅里叶变换DFT(快速傅里叶变换FFT)进行频谱分析》由会员分享，可在线阅读，更多相关《Matlab如何进行离散傅里叶变换DFT(快速傅里叶变换FFT)进行频谱分析（15页珍藏版）》请在金锄头文库上搜索。
1、Matlab 如何进行如何进行利用利用离散傅里叶变换离散傅里叶变换 DFT (快速傅里叶变换快速傅里叶变换 FFT)进行进行频谱分析频谱分析 信号在频域能够呈现出时域不易发现的性质和规律， 傅里叶变换是将信号从时域变换到频域，便于在频域对信号的特性进行分析。离散傅里叶变换 (DFT)，是傅里叶变换在时域和频域上的离散呈现形式， 通俗的说就是将经过采样的有限长度时域离散采样序列变换为等长度的频域离散采样序列， 通过对变换得到的频域采样序列进行适当的换算和处理，可以得到信号的频谱（频率-幅值曲线和频率-相位曲线） 。 离散傅里叶变换 (DFT)的定义为： 210( )DFT ( )( ),01Nj
2、nkNnX kx nx n ekN 式中，x(n)为时域离散采样序列（通常为实数序列） ，N 为时域离散采样序列 x(n)的长度， X(k)为频域离散采样序列 （通常为复数序列） 。 快速傅里叶变换（FFT）是离散傅里叶变换（DFT）的一种快速算法，FFT 的计算结果与 DFT 完全相同，但 FFT 相对于 DFT 减小了计算量、节约计算资源消耗，能够适应在线计算，因此实际 DFT 都是通过 FFT 算法来求得结果。 Matlab 软件自带 fft 函数实现快速傅里变换算法，但是光使用 fft并不能直接得到信号的频谱，还需要解决以下问题： |X(k)|X(k)kk幅值对称幅值对称对应频率对应频
3、率fs/2直流信号幅值直流信号幅值( (对应频率对应频率0Hz)有效频率范围有效频率范围0 0 fs/2相位相反相位相反幅值大小需要幅值大小需要乘以乘以2/N DFT 变换后的 X(k)复数序列幅值、相位图 1) 幅值变换幅值变换：X(k)序列的幅值大小与参与变换的时域序列 x(n)长度 N 有关，变换后的幅值|X(k)|需要乘以 2/N 得到真实幅值； 2) 有效频率区域有效频率区域：X(k)序列由两部分共轭复数序列组成（复数共轭表示幅值相等、相位相反） ，相当于只有一半的复数序列是独立有效的，这部分复数序列对应 0fs/2 的频率区域（fs为时域离散采样序列 x(n)的采样频率） 。 3)
4、 直流信号的处理直流信号的处理：直流信号幅值（对应频率 0Hz）为两部分共轭复数序列在频率 0Hz 处的加和，其真实幅值再乘以 2/N 后还需要再除以 2 得到真实的直流信号幅值。 初学的朋友若不理解上述变换和处理技巧， 很难得到正确的频谱图。为此作者在 fft 函数的基础上，使用 Matlab 开发了函数 DFT.m，通过函数来实现上述幅值变换、有效频率区域和直流信号的处理，能够直接分析出给定离散信号 x(n)的幅值谱和相位谱， 函数简单、 易用、通用性好。 function f,X_m,X_phi = DFT(xn,ts,N,drawflag) % f,X_m,X_phi = DFT(xn
5、,ts,N,drawflag) 离散序列的快速傅里叶变换，时域转换为频域 % 输入 xn为离散序列 为向量 % ts为序列的采样时间/s % N为FFT变换的点数，默认为xn的长度 % drawflag为绘图标识位，取0时不绘图，其余非0值时绘图，默认为绘图 % 输出 f为频率向量 % X_m为幅值向量 % X_phi为相位向量，单位为 % 注意计算出来的0频分量(直流分量应该除以2) 直流分量的符号应结合相位图来确定 % By ZFSwust 2020 % 获取更多Matlab/Simulink原创资料和程序，清关注微信公众号：Matlab Fans 下面结合实例进行演示和分析。 例例1：单
6、频正弦信号单频正弦信号（整数周期采样）（整数周期采样） % Eg 1 单频正弦信号 ts = 0.01; t = 0:ts:1; A = 1.5; % 幅值 f = 2; % 频率 w = 2*pi*f; % 角频率 phi = pi/3; % 初始相位 x = A*cos(w*t+phi); % 时域信号 figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,ts); 结果结果： 结果： 正弦信号频率为2Hz，频谱分析频率为1.98Hz 正弦信号幅值为1.5，频谱分析幅
7、值为1.495 正弦信号相位为60，频谱分析相位为63.32 例例 2：单频正弦信号（单频正弦信号（非非整数周期采样）整数周期采样） % Eg 2 单频正弦信号(非整数周期采样) ts = 0.01; t = 0:ts:1; A = 1.5; % 幅值 f = 1.5; % 频率 w = 2*pi*f; % 角频率 phi = pi/3; % 初始相位 x = A*cos(w*t+phi); % 时域信号 figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,ts); 结
8、果结果： 正弦信号频率为1.5Hz，频谱分析频率为0.99Hz、1.98Hz 正弦信号幅值为1.5，频谱分析幅值为1.034、0.923 正弦信号相位为60，频谱分析相位为160.93、-33.76 总结：总结： DFT 变换后频率序列的最小单位刻度为 fs/N （此例为 1Hz） ，非整数周期采样时关心信号的频率（此例为 1.5Hz）不是频率分辨率fs/N 的正整数倍，那这个频率成分信号会由前后两个正整数倍的频率成分信号（此例为 1Hz 和 2Hz）的线性组合来替代，这就是频谱泄漏现象， 非周期采样时某频率成分信号向两侧频率分辨率正整数倍的频点泄漏。实际频谱分析时并不清楚所关心的频率点精确值
9、，避免此问题的一个解决方法是，取更多的点参加 DFT，即时域序列 x(n)长度 N值取长一些，让频率分辨率 fs/N 很小，以减小频谱泄漏现象。 % Eg 2 单频正弦信号(非整数周期采样) ts = 0.01; t = 0:ts:1; A = 1.5; % 幅值 f = 1.5; % 频率 w = 2*pi*f; % 角频率 phi = pi/3; % 初始相位 x = A*cos(w*t+phi); % 时域信号 figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,
10、ts); 结果结果：频谱泄漏情况大为改善，采样点继续增多时，频谱泄漏会进一步减小。 正弦信号频率为1.5Hz，频谱分析频率主要成分为1.46Hz、 正弦信号幅值为1.5，频谱分析频率主要成分对应幅值为1.41 正弦信号相位为60，频谱分析频率主要成分对应相位为89.5 例例 3：含有直流分量的单频正弦信号含有直流分量的单频正弦信号 % Eg 3 含有直流分量的单频正弦信号 ts = 0.01; t = 0:ts:1; A = 1.5; % 幅值 f = 5; % 频率 w = 2*pi*f; % 角频率 phi = pi/6; % 初始相位 x = 0.5 + A*cos(w*t+phi);
11、% 时域信号,带有直流偏移0.5 figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,ts); 结果结果： 正弦信号频率为5Hz，频谱分析频率为4.95Hz 正弦信号幅值为1.5，频谱分析幅值为1.498 正弦信号相位为30，频谱分析相位为38.66 正弦信号直流分量0.5，频谱分析直流分量为0.51 例例 4：正弦复合信号正弦复合信号 % Eg 4 正弦复合信号 ts = 0.01; t = 0:ts:2; A = 1.5 1 0.5 0.2; % 幅值 f = 3
12、6 9 15; % 频率 w = 2*pi*f; % 角频率 phi = (1:4)*pi/4; % 初始相位 x = -0.5 + A(1)*cos(w(1)*t+phi(1) + A(2)*cos(w(2)*t+phi(2) + A(3)*cos(w(3)*t+phi(3) + A(4)*cos(w(4)*t+phi(4); % 时域信号 figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,ts); 结果结果： 正弦信号频率为3、6、9、15Hz，频谱分析频率为2.
13、985、5.97、8.96、14.93Hz 正弦信号幅值为1.5、 1、 0.5、 0.2， 频谱分析幅值为1.499、 0.989、 0.485、 0.192 正弦信号相位为45、90、135、180，频谱分析相位为50.6、101.5、152.9、210 正弦信号直流分量-0.5，频谱分析直流分量为-0.497 注意注意：频率为 0Hz 时对应的直流信号的幅值的正负号，是通过零频相位来确定的，相位为 0表示幅值为正，相位为 180表示幅值为负。 例例 5：含有随机干扰的正弦信号含有随机干扰的正弦信号 % Eg 5 含有随机干扰的正弦信号 ts = 0.01; t = 0:ts:2; A =
14、 1 0.5; % 幅值 f = 3 10; % 频率 w = 2*pi*f; % 角频率 phi = (1:2)*pi/3; % 初始相位 x = A(1)*cos(w(1)*t+phi(1) + A(2)*cos(w(2)*t+phi(2) + 0.8*(rand(size(t)-0.5); % 时域信号 figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,ts); 结果结果： 正弦信号频率为3、10Hz，频谱分析频率为2.985、9.95Hz 正弦信号幅值为1、0
15、.5，频谱分析幅值为0.978、0.456 正弦信号相位为60、135，频谱分析相位为65.1、139.8 例例 6：实际案例实际案例 load data ts = 0.001; x = Jsd; t = 0:length(x)-1*ts; figure plot(t,x) xlabel(时间/s) ylabel(时域信号x(t) % DFT变换将时域转换到频域,并绘制频谱图 f,X_m,X_phi = DFT(x,ts); 结果结果： 频谱分析主要频率成分为18.996、37.992Hz 频谱分析主要频率成分对应幅值为1.741、1.117 该项目为作者在强振环境下测得加速度信号， 加速度是
16、机械结构周期运动激励产生，需要通过频谱分析获取机械结构周期运动的频率。由于噪声幅度远大于有效信号幅度，信号的信噪比很低，从时域上很难辨别机械结构周期运动的频率。但经过 DFT 后，从频域上可以看出信号的主要频率成分为 19Hz 和其倍频 38Hz， 可以判断机械结构周期运动的频率为 19Hz，38Hz 为结构响应的非线性特性所产生的倍频。 工程上我们还会遇到这样的问题：获取了信号的频谱，希望从信号的频谱来恢复时域信号。这就涉及离散傅里叶逆变换问题，下一篇详谈。 有 Matlab/Simulink 方面的技术问题， 欢迎发送邮件至 讨论。 添加 QQ:944077462，免费获取源程序。 更多 Matlab/Simulink 原创资料，欢迎关注微信公众号：Matlab Fans
