信号处理之快速傅里叶变换(FFT)-通俗易懂
信号处理之快速傅里叶变换FFT-通俗易懂
-
- 历史回顾
- 欧拉方程
(FS) 傅立业级数
(FT) 傅立业转换
(DFS) 离散法立业级數
(DTFT) 離散時间法立業轉換
(DFT) 離散法立業轉換
(FFT) 快速法立業轉換（FFTransform）
(MATLAB中常用的工具) MATLAB常見的FFTransform工具包
(在FFT中的常见问题分析) FFT中的常見問題探討
- 欧拉方程
- 历史回顾
历史溯源
相信许多人已经了解了傅里叶变换这一概念。然而，在理解上仍存在模糊不清的状态。今天作者将围绕这一重要数学工具的发展历程及其应用意义展开深入探讨。
约瑟夫·傅里叶（1768年3月21日—1830年5月16日），作为一位多才多艺的数学家兼物理学家，在研究中首次提出了著名的傅里叶级数，并将其成功应用于热传导理论和振动理论分析中。
欧拉公式
欧拉公式属于复分析领域的著名公式，在将三角函数与复指数函数连接方面具有重要意义，并因由莱昂哈德·欧拉命名而被公认为最重要的十个数学公式的之一。欧拉公式已被选为世界上最伟大的十个公式之一。对于任意实数值\phi而言：
e^{i\phi} = \cos(\phi) + i\sin(\phi) 其中 e 代表自然对数的底数；i 是虚数单位；而 \cos 和 \sin 分别对应余弦和正弦等三角函数；参数 \phi 则以弧度表示的角度值。如图所示，在复平面中；点 e^{i\phi} 可以表示圆上任意一点；通过其实部与虚部的组合，则可与该平面中的任一点一一对应起来。因此；基于欧拉公式的指数形式即可构建整个复平面体系；此外；傅里叶级数经频域转换后呈现的就是一种基于复数形式的表现形式；因此；通过一个指数形式就能实现两者之间的统一。
傅里叶级数(FS)
（数学定义） 傅立叶级展开数：它指出任何周期连续的函数 都可以用正弦函数和余弦函数构成的无穷级数来表示。这种级数之所以被称为傅立叶级数，是因为它是由傅里叶提出的。傅里叶选择正弦函数与余弦函数作为基函数，原因在于它们是正交的。
在代数中，我们用基底的线性组合可以表示任意的一组向量；在函数中，我们也可以用基底表示任意的函数。
在代数中，我们希望基底是正交的，方便我们寻找线性组合的系数，例如正交单位向量a = [1,0,0], b = [0,1,0], c = [0,0,1], 那么对于空间中中任意一个向量d = [4,2,1] = 4 a + 2 b + c 来表示；在函数分解中，我们也希望基底是正交的。傅里叶用三角函数做了基底，比如{sin t, cos t, sin 2t, cos 2t, sin 3t, … sin nt, cos nt}。可以证明这组三角函数基底是正交函数基底， 同理任意一个周期函数可以通过三角函数这组基底表示 例如： f(x) = a*sin t + b * cos t。
那么，何为内积何为正交呢？ 两个平面向量正交的时候时垂直的，写成向量乘法就是\vec{a} \cdot \vec{b} = 0。在学习了线性代数后，我们把它写成了 \bar{X_{1}} * \bar{X_{2}}^{T} = 0。这里的向量可以是任意维数的，比如\left ( X_{1}, X_{2}, X_{3}...X_{n} \right ) 。上面的点乘被称为求取向量的内积，即对应元素的求积累加，正交就是内积为0。 那么，同理，一组正交函数基底的内积的正交性可以用积分形式表示：\int_{a}^{b} f_{1}(x) * f_{2}(x) = 0，两函数正交，可把它叫做内积积分为0。这里函数基底可以是任意维数的，比如：\left ( f_{1}(x) , f_{2}(x), f_{3}(x)... f_{n}(x) \right )。
在代数中，向量在一个基底上的线性组合的系数可以用内积得到；在几何中，向量在一个基底上的线性组合的系数就是向量在该基底的投影，即内积就是投影 。在函数中，我们用内积积分可以求得相应地系数，函数在一个基底上的系数也是在该基底上的投影。
例如： 向量 \vec{a} 在 一组基底 \left ( X_{1}, X_{2}, X_{3}...X_{n} \right ) 中的系数分别是 \vec{a} * X_{1}， \vec{a} * X_{2}，\vec{a} * X_{3}… \vec{a} * X_{n}， 其中\vec{a} * X_{1} = \left | \vec{a} \right | * \left | X_{1}\right | * \cos (\Theta )，这就是投影。
同理，傅里叶级数展开的系数 ：对于任意函数f(x), 在正弦基底\left ( \sin(\frac{2\pi }{T}), \sin(\frac{4\pi }{T}), \sin(\frac{6\pi }{T})...\sin(\frac{2n\pi }{T}) \right )和余弦基底\left ( \cos (\frac{2\pi }{T}), \cos (\frac{4\pi }{T}), \cos (\frac{6\pi }{T})...\cos (\frac{2n\pi }{T}) \right )中展开得到，其中\omega0 = \frac{2\pi }{T},即我们所说的角频率或频率或最小分辨率，\frac{2n\pi }{T} 为n倍倍频率。对于任意一个三角函数：A *\sin (\omega0 *x + b)，A即为振幅，\frac{2\pi }{T} = \omega0是角频率，b为初相，\omega0 *x + b为相位，这样，傅里叶级数展开后，幅频特性和相频特性就可以计算出来了。如下公式，可试着将n=1,2,3,…,N代入公式计算推演一下，由此可以看出傅里叶级数展开的结果是非周期离散的 。
用欧拉公式换算成指数形式：
需要注意的是，在公式\frac{2\pi nx}{T} = \omega_0 \times nx中代表的是一个相位值。其中\omega_0代表角频率或者说是最小的频率分辨率参数。而\omega_0 \times n则表示该基频值的一个整数倍频数关系。
傅里叶变换(FT)
在上一部分傅里叶级数(FS)中,我们是通过三角函数或指数形式来表示周期连续函数的方法,但对于非周期连续函数该如何处理呢,而且非周期连续函数比周期连续函数更为普遍,这也是我们引入傅里叶变换(FT)的原因所在。所谓傅里叶变换就是将傅里叶级数的概念推广到非周期连续函数领域中去。
为了实现这一目标,我们借助极限思想进行处理:假设非周期连续函数在整个定义域内都包含在一个无限大的周期区间内,这样该区间就不再具有任何重复性特征了,从而我们可以将其视为一个具有无穷大周期的周期连续函数来进行分析研究。这种情况下,我们可以通过以下数学表达式来进行形式上的描述:当连续非周期函数的周期T趋向于无穷大时,f(x)就可以被看作是一个以无穷大为周期的周期连续函数了。
具体来说,对于这样的情况,我们可以得到如下的数学表达式:当T→∞时,Cn=1/T∫x0x0+Tf(x)e{-i2πTx}dx=1/T∫0Tf(x)e{-inω0x}dx=C(nω0),其中ω0=2π/T随着T趋向于无穷大而趋向于零值,n是一个整数变量趋于无穷大的过程则会导致nω0这个离散角频率变得不再离散而是变成了一个连续的角频率变量并且其取值范围是有限大的值为此我们将nω0定义为新的变量ω即令nω=ω这样原来的积分区间就被扩展到了整个实数轴上即有如下的积分表达式
C(ω)=∫−∞+∞f(x)e^{-iωx}dx
同样地,f(x)也可以表示为f(x)=1/(2π)∫−∞+∞C(ω)e^{iωx}dx
通过以上推导可以看出,C(ω)实际上是以频率为自变量的一个复变函数,而其结果也自然地成为了非周期连续频谱的表现形式
离散傅里叶级数(DFS)
考虑一个 period 为 T_0 的 continuous periodic signal f(x), 在该 period 内均匀采样 N 个 point, 即 x_0, x_1, ..., x_{N-1}, 其中 x_m = f(m*T) (m=0,...,N-1), 这里的 T 表示采样间隔时间（即 T = T_0/N）。
由此可得，在时域上表现为 discrete 且 periodic 的序列 x[n]（其中 n=..., -2, -1, 0, 1, ..., N-1）满足关系式 x[n + N] = x[n]（对于任意整数 n）。
离散时间傅里叶变换(DTFT)
上一部分我们掌握了离散周期信号的处理方法接下来我们探讨离散非周期信号的处理方案DTFT方法通过均匀采样将连续时间非周期信号转化为离散序列进而求取其傅里叶变换从而得到离散非周期信号的频域特征这一过程也被称为离散时间傅里叶变换从理论上讲通过傅里叶变换理论可以得到：将函数f(x)在区间(-∞ +∞)内进行积分运算即可获得其频域表示C(ω) = ∫{−∞}^{+∞} f(x)e^{-iωx} dx经过采样后我们有x[n] = f(nT)其中T为采样间隔N为采样点数这样就可以将上述积分表达式转化为求和形式即C(ω) = Σ{n=-∞}^{+∞} x[n]e^{-iωnT}这表明对连续信号进行均匀采样后其频谱特性将会发生特定的变化
离散傅里叶变换(DFT)
因为对离散非周期信号进行变换后得到的是连续周期性的频谱，在现实中我们无法直接应用这种频谱。因此我们需要将这种连续周期的频谱在频域上进行采样处理使其变为离散化的形式从而引出了离散傅里叶变换（DFT）。在此过程中我们将频率ω被以2π/N为间隔进行离散化处理进而使得DFT得以实现其数学表达式为C(ω)=∑{-∞}{+∞}x[n]e{-iωn}通过这一过程我们能够获得所需的数字频谱信息其中C[k]表示第k个频率点上的复数幅度值其计算公式为C[k]=∑{n=0}{N-1}x[n]e{-i(2πk/N)n}其中k=0,1,2,…,N-1分别对应于不同频率点上的采样值。同样地为了恢复原始信号我们也需要对这些复数幅度值进行反变换其计算公式为x[n]=∑_{k=0}{N-1}C[k]e{i(2πk/N)n}通过这样的正反变换过程我们可以完全恢复出原始信号的时间域信息并从中提取出幅值信息以及相位信息这些信息共同构成了完整的频谱分析结果如图所示展示了不同频率点上的幅值分布情况。
总结：
- 各种信号时域和频域的关系(连续对应非周期，离散对应周期)
| 时域 | 频域 |
|---|---|
| 连续 、非周期性 | 非周期性、连续 |
| 连续 、周期性 | 非周期性、离散 |
| 离散 、非周期性 | 周期性、连续 |
| 离散 、周期性 | 周期性、离散 |
各种信号的图形实例
快速傅里叶变换(FFT)
在之前的章节中，我们深入学习了离散傅里叶变换（DFT）的相关知识。然而，在公式推导中
基于公式的结果显示, 每个X[k]的计算涉及N次复数乘法操作和N-1次复数加法操作, 而每个X[k]点所需的总运算量则分别为4N次实数乘法操作以及2N+2(N-1)次实数加法操作. 这些运算量随着数据点数量的增长呈平方关系上升, 因此, 为了提高运算效率, 必须采用高效的快速傅里叶变换算法（FFT）. 其中基2 FFT 算法因其在单片机平台上的广泛应用而成为主流选择. 该算法的基本思想是将输入信号按照时间顺序进行偶奇分解, 并通过递归迭代的方式逐步降维至二维甚至一维数据进行处理. 这种分层迭代的过程可以用"按时间抽取"来形象地描述: 对于一个长度为N=2^{M}的数据序列进行DFT变换, 首先将其分解为两个长度为N/2的小序列分别进行DFT, 再进一步分解直至只剩下一个数据点. 最终通过蝶形运算组合得到最终结果. 如果原始信号长度不是严格的2的幂次数, 可以在尾部补零使其满足条件.
4点基2FFT实例计算流程：
8点实例基2FFT计算流程：
MATLAB中常用的FFT工具
其中 fft(X), fft(X,N) 中,X代表输入信号序列. 在MATLAB中, fft函数基于DFT算法进行运算. 即, fft(X) 默认以X信号序列长度作为计算点数. 同样地, 可以通过 fft(X,N) 指定所需的FFT点数.
close all;
clear ;
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*100*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,100,200的叠加信号
N = length(x);
% FFT变换
Y = fft(x);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y).^2/N;
% 绘制频谱
figure
plot(f, power(1:N/2+1));
title('FFT of Windowed Signal');
xlabel('Frequency (Hz)');
ylabel('power');
AI写代码
[pxx,f] = pwelch(X,window,noverlap,NFFT,fs)主要用于分段周期图法来估算功率谱密度（PSD）。
其中：
- x是一维的时间序列数据；
- window用于表示每个周期图窗口所包含的时间序列数据点数量；
- 选择合适的窗宽建议参考相关公式；
- n_overlap参数表示相邻两个周期图之间的重叠点数；
- N_FFT参数设置应大于每一段的数据点数，并且最好选择比其大的最小二进制数值以提高频域分辨率；
- 采样频率fs决定了频域分析的最大频率值（f_max = fs/2）；
- 功率谱密度（PSD）是一种用于分析信号频域特性的工具，在工程实践中具有重要应用价值；
PSD通过将时域信号转换至频域来揭示信号能量随频率分布的情况；通俗而言，在不同频率区间内信号的变化程度可以通过PSD进行量化评估；这对于深入理解信号特性及提取关键信息非常有帮助。
close all;
clear ;
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*100*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,100,200的叠加信号
N = length(x);
K = 2; %把1000个点分为2段
M = N/K; %M为每一段加窗的点数
[pxx,fxx] = pwelch(x,hanning(M),0,N,Fs);
plot(fxx,pxx);
grid;
xlabel('Frequency ,Hz')
ylabel('pwelch函数估计的PSD ,dB')
AI写代码
3. [pxx,f] = plomb(x,t)。其中x代表输入的信号序列,t则对应时间序列变量。
隆布-斯卡格尔方法在汽车工业、通信系统以及医学研究等领域的应用尤为广泛。
该方法特别适用于处理非均匀采样数据，在现代信号处理中具有重要价值。
对于均匀采样的信号处理而言,传统频谱分析方法确实是一种标准工具。
然而,当遇到非均匀采样情况时,直接应用这类方法将导致结果失真或计算误差。
为了克服这一挑战,隆布-斯卡格尔算法提供了一种直接处理非均匀采样数据的有效途径。
该算法已在plomb函数模块中得到实现,但其内部工作原理仍需进一步深入研究。
在天文学研究领域,这一方法的应用前景非常广阔。
FFT中常见的问题
- 频率混叠
对连续信号进行等时间采样时，如果采样频率不满足采样定理，采样后的信号频率就会发生混叠，即高于奈奎斯特频率（采样频率的一半）的频率成分将被重构成低于奈奎斯特频率的信号。这种频谱的重叠导致的失真称为混叠，也就是高频信号被混叠成了低频信号。采样定理 ，又称香农采样定理，奈奎斯特采样定理，只要采样频率大于或等于有效信号最高频率的两倍，采样值就可以包含原始信号的所有信息，被采样的信号就可以不失真地还原成原始信号。如下图红色是真实信号，蓝色为采样信号，由于不满足采样定理，严重失真。
为了防止频谱混叠现象的发生，在对原始信号进行降采样操作的过程中，则必须先通过低通滤波将信号调制至有效频段范围（即小于等于奈奎斯特频率），之后再按规定的抽样间隔进行抽样以获得目标频段下的数字序列。
close all;
clear ;
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*30*t) +0.9*sin(2*pi*100*t)+ 0.3*sin(2*pi*200*t); % 生成频率为30,100,200的叠加信号
N = length(x);
% FFT变换
Y = fft(x);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y).^2/N;
figure
subplot(321)
plot(x)
title('1000Hz采样信号');
xlabel('时间');
ylabel('幅值');
subplot(322)
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('功率');
fsr = 100;
x1 = resample(x, fsr,Fs);
t1 = 0:1/fsr:1-1/fsr;
N = length(x1);
% FFT变换
Y = fft(x1);
% 计算频率轴
f = fsr*(0:(N/2))/N;
%功率
power = abs(Y).^2/N;
subplot(323)
plot(x1)
title('直接将采样到100Hz信号');
xlabel('时间');
ylabel('幅值');
subplot(324)
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('功率');
[bb, aa] = butter(4, 50*2/Fs,'low'); % 50Hz低通滤波
data1 = filter(bb, aa, x);
fsr = 100;
x1 = resample(data1, fsr,Fs);
t1 = 0:1/fsr:1-1/fsr;
N = length(x1);
% FFT变换
Y = fft(x1);
% 计算频率轴
f = fsr*(0:(N/2))/N;
%功率
power = abs(Y).^2/N;
subplot(325)
plot(x1)
title('先50Hz低通滤波再降采样到100Hz信号');
xlabel('时间');
ylabel('幅值');
subplot(326)
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('功率');
AI写代码- 频谱泄露
所谓频谱泄露是指信号频谱中各频率之间相互影响导致测量结果偏离实际值的同时在各频率邻近两侧其他频率点上出现一些幅值较小的假谱其主要原因有两个第一是由于最小分辨频率不是Fs/N的整数倍其中Fs为采样率N为采样点数目因此离散傅里叶变换DFT仅能输出Fs/N整数倍频率上的功率当输入信号不在这些特定频率点时DFT输出将无法准确反映相应频率信息从而导致能量分布不准确具体的泄漏分布情况取决于所采用窗函数在连续域中的傅里叶变换特性第二个原因是施加了矩形窗如果直接对信号进行快速傅里叶变换FFT等价于对该信号施加了一个矩形窗这会导致能量向非预期的旁瓣频率扩散从而加剧频谱泄露问题为了减轻频谱泄露的影响可以通过增大窗函数宽度或者采用加权型窗函数例如汉明窗汉宁窗和高斯窗等这些加权函数能够有效减少旁瓣的能量含量并改善频谱分析结果下图展示了对包含30Hz 100Hz和200Hz叠加信号进行分析的情况原始采样率为1000Hz末尾补零前直接进行FFT分析会观察到明显的旁瓣现象而如果在补零之前先施加一个合适的窗函数则可以使旁瓣能量得到显著抑制从而提高分析精度因此在执行补零操作之前必须先选择一个合适的窗函数进行处理
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*100*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,100,200的叠加信号
N = length(x);
% FFT变换
Y = fft(x,N);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y)/N;
figure
subplot(321)
plot(x);
title('1000Hz信号');
xlabel('时间');
ylabel('幅度');
subplot(322);
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('幅值');
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*100*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,100,200的叠加信号
x = [x,zeros(1,length(t))];
N = length(x);
% FFT变换
Y = fft(x,N);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y)/N;
subplot(323)
plot(x);
title('1000Hz信号末尾补零');
xlabel('时间');
ylabel('幅度');
subplot(324);
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('幅值');
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*100*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,100,200的叠加信号
N = length(x);
% % 加汉明窗
window = hanning(N);
x = x .* window';
x = [x,zeros(1,N)]; %末尾补length(x)个零
N = length(x);
% FFT变换
Y = fft(x,N);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y)/N;
subplot(325)
plot(x);
title('1000Hz信号先加窗后末尾补零');
xlabel('时间');
ylabel('幅度');
subplot(326);
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('幅值');
AI写代码- 栅栏效应
在进行DFT的过程中，最后需要对信号的频谱进行采样。经过这种采样所显示出来的频谱仅在各采样点上，而不在此类点上的频谱都显示不出来，即在其他点上有重要的峰值也会被忽略，这就是栅栏效应。也就是说，由于对频谱离散采样，有些重要的频率点未被采样到，结果丢失了重要的频率信息。就像栅栏一样，栅栏缝隙越宽，能拦住的东西越少，栅栏缝隙越窄，能拦住的东西越多，所以要想化解这种问题就要提高频率分辨率。即可通过信号加窗末尾补零的方法提高采样点数，进而提高频率分辨率。
不管是时域采样还是频域采样，都有相应的栅栏效应。只是当时域采样满足采样定理时，栅栏效应不会有什么影响。而频域采样的栅栏效应则影响很大，“挡住”或丢失的频率成分有可能是重要的或具有特征的成分，使信号处理失去意义。减小栅栏效应可用提高采样间隔也就是频率分辨力的方法来解决。间隔小，频率分辨力高，被“挡住”或丢失的频率成分就会越少。但会增加采样点数，使计算工作量增加。
举例：如下图所示采样率为1000，采样了1000个点，生成频率为50,70.5,200的叠加信号，由于第一行信号的频率分辨率为Fs/N=1Hz，因此对于在频率为70.5Hz的点在FFT结果后就会出现畸变，丢失了频率成分，发生了栅栏效应。第二行，通过末尾补零，有效缓解了栅栏效应，但有些频率泄露。第三行是先加hamming窗再补零，效果明显改善了很多。
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*70.5*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,70.5,200的叠加信号
N = length(x);
% FFT变换
Y = fft(x,N);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y)/N;
figure
subplot(321)
plot(x);
title('1000Hz信号');
xlabel('时间');
ylabel('幅度');
subplot(322);
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('幅值');
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*70.5*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,70.5,200的叠加信号
x = [x,zeros(1,length(t))];
N = length(x);
% FFT变换
Y = fft(x,N);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y)/N;
subplot(323)
plot(x);
title('1000Hz信号末尾补零');
xlabel('时间');
ylabel('幅度');
subplot(324);
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('幅值');
Fs = 1000; % 采样频率为1000
t = 0:1/Fs:1-1/Fs; % 时间向量
x = 0.7*sin(2*pi*50*t) +0.9*sin(2*pi*70.5*t)+ 0.3*sin(2*pi*200*t); % 生成频率为50,70.5,200的叠加信号
N = length(x);
% % 加汉明窗
window = hanning(N);
x = x .* window';
x = [x,zeros(1,N)]; %末尾补length(x)个零
N = length(x);
% FFT变换
Y = fft(x,N);
% 计算频率轴
f = Fs*(0:(N/2))/N;
%功率
power = abs(Y)/N;
subplot(325)
plot(x);
title('1000Hz信号先加窗后末尾补零');
xlabel('时间');
ylabel('幅度');
subplot(326);
plot(f, power(1:N/2+1));
title('FFT结果');
xlabel('频率');
ylabel('幅值');
AI写代码旁瓣效应
在时域中可以选择不同的窗函数用于截断信号，在这种情况下可以使相邻频点之间的旁瓣高度趋向于零的同时将能量主要集中在主瓣区域。然而加窗操作本身会引入误差，并且不同类型的窗函数会导致能量分配上的差异；因此选择合适的窗函数取决于输入信号的特性以及测试人员所关注的问题领域（如低频或高频分析）。我们通常面临一个权衡问题：矩形窗具有较窄的主瓣宽度（频率分辨率较高）但其旁瓣衰减能力较弱；而高斯型窗等其他类型的窗函数虽然能有效抑制旁瓣振幅（降低频谱拖尾）但在某些情况下可能会牺牲主瓣宽度以换取较低的频率分辨率（如图所示）。对于非常邻近的频率分量分析（例如频率分辨率为0.1Hz），若选择宽主瓣的矩形窗可能导致1Hz与1.1Hz等频点的表现不够明显；此时建议采用主瓣宽度较窄的其他类型 window 进行分析以提高区分度
总结：
在采样过程中x(t)时域离散化的过程中可能会出现频谱混叠现象（即f_s/2$以下的高频信号会被折叠到低频区域），特别在进行降采样操作时必须确保不会引入混叠现象（即必须先对信号进行低通滤波器处理再进行降采样）。**为了减少频谱泄露的影响需合理选择适当的窗函数来进行加窗处理；**为了消除栅栏效应应在完成加窗处理后对信号序列补充足够的零点使其达到所需的长度要求这也是一种常见的处理手段。