本文纯原创，有不对的地方敬请指正，有不懂的地方可以讨论，搬运请注明来源，谢谢。
目录
一、问题
在研究matlab的FFT函数的时候发现了如下问题：对于信号
y
=
e
j
2
π
f
1
t
+
e
j
2
π
f
2
t
+
e
j
2
π
f
3
t
y=e^{j2\pi f_1t}+e^{j2\pi f_2t}+e^{j2\pi f_3t}
y=ej2πf1t+ej2πf2t+ej2πf3t
其中
f
1
=
500
H
z
f_1=500Hz
f1=500Hz，
f
2
=
505
H
z
f_2=505Hz
f2=505Hz，
f
2
=
1010
H
z
f_2=1010Hz
f2=1010Hz。当我用
F
s
=
5120
H
z
F_s=5120Hz
Fs=5120Hz的采样率对其进行采样得到512个离散数据点，画出其频谱为：
问题一：频谱本应该是在500Hz、505Hz和1010Hz三处有峰，但是频谱上却没有显示出505Hz的信号。
问题二：由时域表达式可知，三个不同频点的信号功率一致，都为1。但是幅度谱上却无体现。
二、频谱分辨率、栅栏效应及频谱泄露
在简介栅栏效应之前，我们要分清DTFT、DFT和FFT三者的关系，可以参考我的这篇文章。
我们不妨细看DFT这个步骤，我们在matlab里调用FFT函数就是采用的DFT算法（当FFT点数是2的整数次幂的时候，就可以用快速傅里叶变换）。
2.1频谱分辨率
由上图可知，我们最终得到的时域信号是以采样周期 T = 1 F s T=\frac{1}{Fs} T=Fs1为间隔，频谱中的谱线以 2 π N \frac{2\pi}{N} N2π间隔，转换到数字上则有： Δ f = F s N \Delta f=\frac{Fs}{N} Δf=NFs，在本例中频谱分辨率则是 Δ f = 5120 512 = 10 H z \Delta f=\frac{5120}{512}=10Hz Δf=5125120=10Hz。
2.2栅栏效应
我们可以想象，间隔10Hz的谱线，应该是会跳过我们的505Hz这个频点的信号。就好像通过栅栏观察频谱一样，因此叫做栅栏效应。
那么难道我们就观察不到这个频点的信号了吗？如果我们漏掉了505Hz这个频点的信号，那我们500Hz和1010Hz两个频点的信号功率不应该是一致的吗，怎么会幅度不一样呢？
2.3频谱泄露
由上图可知，DFT的最后一个步骤是加窗，将一个无限长的信号截断为有限长的信号。在频域上，就是卷积一个Sa函数。
Sa函数的特性是能量集中在主瓣，少部分能量分布在旁瓣。这样原本各个频点的冲激信号经过卷积变为了Sa函数，能量也从集中分布在某个频点变为大部分集中，小部分泄露的情况。这就是我们说的频谱泄露。
回归到这个例子中，505Hz频点处的信号能量理应集中分布在505Hz附近的主瓣中，少部分泄露到其他地方去了。
可是505Hz谱线最后由于频谱分辨率的原因并没有被采集到，而500Hz的谱线不仅集中了500Hz信号的能量，还因为离505Hz太近，包含在其主瓣内，分到了505Hz频点信号的能量，因此幅度变得特别大。而1010Hz离得很远，只分到了一点点的能量，所以它的幅度也不是整数。（具体该是多少我们下一章节讨论）
由此可见，如果需要缓解频谱泄露的问题，如何选择一个好的窗函数是很重要的。
2.4如何提高频谱分辨率
如要提高频谱分辨率，还是从公式考虑
Δ
f
=
F
s
N
\Delta f=\frac{Fs}{N}
Δf=NFs
要缩小
Δ
f
\Delta f
Δf，采样率
F
s
Fs
Fs一般是不做改动，而是增大
N
N
N。通常大家会提到**“补零”、“重复数据”和“增加采样点”**三种方法，我们一一讨论
2.4.1补零是否可以提高分辨率？
频谱分辨率公式里的
N
N
N究竟是什么呢？是采样点个数还是FFT点数？首先要分辨频谱密度和频谱分辨率两个不同的概念，经常有人对此混淆不清。
频谱密度=
F
s
N
f
f
t
\frac{Fs}{N_{fft}}
NfftFs，
N
f
f
t
N_{fft}
Nfft是FFT点数。
频谱分辨率=
F
s
m
i
n
{
N
,
N
f
f
t
}
\frac{Fs}{min\{N,N_{fft}\}}
min{N,Nfft}Fs，
N
N
N是采样点个数。
我们讨论以下三种情况：
采样点个数=FFT点数
易知，当采样点个数与FFT点数相同时，频谱密度=频谱分辨率= F s N f f t \frac{Fs}{N_{fft}} NfftFs= F s N \frac{Fs}{N} NFs。
采样点个数>FFT点数
当采样的数据点过多时候，如果做大点数的FFT变换，运算量会很大。一般会选择1024点，16384点这样的点数进行快速傅里叶变换以提高计算效率。
此时频谱密度=频谱分辨率=
F
s
N
f
f
t
\frac{Fs}{N_{fft}}
NfftFs。但是增加FFT点数可以继续增大频谱分辨率，直到
F
s
N
\frac{Fs}{N}
NFs为止，这也可以说是频谱分辨率的上限。
采样点个数< FFT点数
这就是我们常说的补0的情况。在这个例子中，我们用
F
s
=
5120
H
z
F_s=5120Hz
Fs=5120Hz的采样率对其进行采样得到512个离散数据点，结果没有分出500Hz和505Hz两个频点的信号。那我们在512个数据点后面补512个0，然后做1024点的FFT，这样频谱分辨率会不会变成5Hz，就可以区分出两个相邻频点的信号呢？
由上图可知，我们仍然无法区分出500Hz和505Hz两个相邻频点的信号，但是我们的谱线明显更平滑了。
其实，从数学的角度出发，我们可以知道，时域补0对应的是频域插值，插值只会让谱线更平滑，而不会引入新的频率分量。从公式的角度出发，增大了FFT的点数，频谱密度增大了，而频谱分辨率没有提升，因此还是无法区分相邻频点。
如果了解DFT和DTFT的关系就更容易理解， 补零越多，就是在DTFT连续的频谱上采样得到更密集的谱线而已。看我上面那个链接，再结合这篇文章可以很透彻的理解
但是补0也不是毫无作用的，此时看频谱的栅栏更密集了，更容易进行观察和分析。通常我们认为补0是无害的，可以适当补0以完成快速傅里叶变换。例如把500个点补到512个点进行FFT。
（题外话，时域插值是否可行）
既然时域补零对应频域插值，那时域插值对应的就是频域的补0，自然也无法引入新的频域分量。这是上采样的工作原理，有兴趣的可以看我这篇文章
2.4.2重复数据是否可以提高频谱分辨率
首先声明，本文不讨论周期信号进行整周期采样的情况。
一般来说我们采样得到的信号都是非周期的，不管是对无限长非周期信号的采样，还是对周期信号的截断采样。
回到我们的例子，我们我们用
F
s
=
5120
H
z
F_s=5120Hz
Fs=5120Hz的采样率对其进行采样得到512个离散数据点，然后把它重复一遍，拼凑得到1024个数据点进行FFT
可以看出，确实多了一个频率分量，区分出来两个相邻频点的信号了。这是因为重复一遍信号，相当于增加了实际有用信息的采样点个数，因此不仅可以缓解频谱泄露的问题，还可以增大频谱分辨率。此时的频谱分辨率 Δ f = 5120 1024 = 5 H z \Delta f=\frac{5120}{1024}=5Hz Δf=10245120=5Hz。
不过这种方法也有缺点，它伪造了一个与实际信号可能并不相符的新的信号出来，这样得到的频谱与我们实际信号的频谱有所区别。只有当我们无法获取更多采样点的时候才考虑这种方法。
2.4.3增加采样点是否可以提高频谱分辨率
先说结论，这是唯一正确的方法。增加采样点个数相当于加窗的时候加了一个更长的窗，不仅有效缓解了频谱泄露的问题，也大大提升了频谱分辨率。只是很多时候受限于存储条件和计算性能的限制，我们无法采样得到足够多的数据点进行FFT变换。
回到我们的例子，我们我们用
F
s
=
5120
H
z
F_s=5120Hz
Fs=5120Hz的采样率对其进行采样得到1024个离散数据点进行FFT，此时的频谱分辨率
Δ
f
=
5120
1024
=
5
H
z
\Delta f=\frac{5120}{1024}=5Hz
Δf=10245120=5Hz。恰好可以分辨出500Hz和505Hz两个相邻频点信号。
而且此时频谱泄露的问题得到了有效缓解，各个频点的信号的功率没有泄漏到其他频点，所有频点信号的功率一致，幅度谱上都显示1024。
由此我们解决了问题一，那么问题二究竟是怎么回事呢？
三、FFT幅度谱数值矫正问题
我们重新来看我们的信号
y
=
e
j
2
π
f
1
t
+
e
j
2
π
f
2
t
+
e
j
2
π
f
3
t
y=e^{j2\pi f_1t}+e^{j2\pi f_2t}+e^{j2\pi f_3t}
y=ej2πf1t+ej2πf2t+ej2πf3t
这是由三个平均功率为1的复单音信号合并而成的一个混合信号。
我们对其采样1024点做1024点FFT，以及采样2048点做2048点FFT得到两个频谱。
我们看幅度，做1024点FFT的三根谱线的幅度都是1024，而做2048点FFT的三根谱线的幅度是2048。由帕萨瓦尔定理我们可知
信号在时域的总能量等于信号在频域的总能量，即信号经傅里叶变换后其总能量保持不变，符合能量守恒定律。
明明他们都是平均功率为1的信号啊，怎么会这样呢？我们还是要从FFT底层入手。
要知道，非周期信号FFT算的是傅里叶变换，周期信号FFT算的是傅里叶级数（因为周期信号傅里叶变换是
δ
\delta
δ冲激函数，没有意义）
回到我最开始引用的那篇文章里的这个图，由于非周期无限长的信号没有办法求得离散频谱，因此我们是对其做一个周期延拓得到周期信号，然后求得离散频谱。
对一个周期延拓的信号，求的是其傅里叶级数。而我们本来要求的是非周期信号的傅里叶变换。那么DFS和DFT求得东西是一样的吗？
我们看公式：
离散周期傅里叶级数DFS
x
[
n
]
=
∑
k
=
0
N
−
1
a
k
e
j
k
(
2
π
/
N
)
n
x[n]=\sum\limits_{k=0}^{N-1}a_k e^{jk(2\pi/N)n}
x[n]=k=0∑N−1akejk(2π/N)n
a
k
=
1
N
∑
n
=
0
N
−
1
x
[
n
]
e
−
j
k
(
2
π
/
N
)
n
a_k = \frac{1}{N} \sum\limits_{n=0}^{N-1} x[n]e^{-jk(2\pi/N)n}
ak=N1n=0∑N−1x[n]e−jk(2π/N)n
离散傅里叶变换DFT
x
[
n
]
=
1
N
∑
k
=
0
N
−
1
X
[
k
]
e
j
k
(
2
π
/
N
)
k
n
x[n]=\frac{1}{N}\sum\limits_{k=0}^{N-1}X[k] e^{jk(2\pi/N)kn}
x[n]=N1k=0∑N−1X[k]ejk(2π/N)kn
X
[
k
]
=
∑
n
=
0
N
−
1
x
[
n
]
e
−
j
k
(
2
π
/
N
)
k
n
X[k] = \sum\limits_{n=0}^{N-1} x[n]e^{-jk(2\pi/N)kn}
X[k]=n=0∑N−1x[n]e−jk(2π/N)kn
仔细比较就可以发现，DFT算出来的X[k]的幅度是 a k a_k ak的N倍，Matlab中FFT使用的就是这个算法，所以实际的幅值大小要乘以1/N。
回到我们的例子，我们在1024点FFT频谱中得到的幅度是1024，乘以1/N得到1，所以各频点信号的平均功率为1，与时域功率一致，符合帕萨瓦尔定理。
但是为什么有的文章说要乘以2/N那？这只是单边谱和双边谱的区别，DFT算出来的是双边谱，换算到单边谱幅度要加倍，注意直流不加倍，仅此而已。
在进行谱分析的过程中我们其实一般不会注意到幅度具体的值，只会关注该点有无信号以及信号强度相对大小，但是若能够了解其幅度变换关系，相信更有利于大家的分析！
扫一扫
569
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
