信号与系统 05：DFT、FFT、窗与频谱标定

完整入口 · 上一章：傅里叶与采样 · 下一章：工程应用

1. 四个名字各解决什么问题

名称	时间输入	频率输出
CTFT	连续时间信号	连续角频率
DTFT	整条离散序列	连续且 $2\pi$ 周期的离散角频率
DFT	有限 $N$ 个点	$N$ 个频率取值
FFT	同 DFT	同 DFT，算法计算更快

本章约定正变换不缩放、逆变换除以变换长度。软件可能采用其他约定，须明确归一化；可对照 NumPy FFT 文档。

2. 从有限序列得到 DFT

把实际记录 $x[0],\ldots,x[N-1]$ 在记录外补零，DTFT 是：

$$X_{\mathrm d}(e^{j\Omega})=\sum_{n=0}^{N-1}x[n]e^{-j\Omega n}.$$

在 $\Omega_k=2\pi k/N$ 读取 $N$ 个点：

$$\boxed{X[k]=\sum_{n=0}^{N-1}x[n]e^{-j2\pi kn/N}}, \quad k=0,\ldots,N-1.$$

这就是 $N$ 点 DFT。它也可看作 $N$ 周期序列傅里叶系数的 $N$ 倍；频域乘法因此对应循环卷积，需要与非周期的线性卷积分开。

3. 为什么逆变换除以 $N$

令 $W_N=e^{-j2\pi/N}$。计算正交和：

$$S_{m-n}=\sum_{k=0}^{N-1}e^{j2\pi k(m-n)/N}.$$

若 $m=n$，每项为 1，和为 $N$。若 $m\ne n$，令 $q=e^{j2\pi(m-n)/N}\ne1$，但 $q^N=1$：

$$S_{m-n}=\frac{1-q^N}{1-q}=0.$$

因此：

$$\begin{aligned} \frac1N\sum_kX[k]e^{j2\pi km/N} &=\frac1N\sum_k\sum_nx[n]e^{j2\pi k(m-n)/N}\\ &=\frac1N\sum_nx[n]S_{m-n}\\ &=x[m]. \end{aligned}$$

所以：

$$\boxed{x[n]=\frac1N\sum_{k=0}^{N-1}X[k]e^{j2\pi kn/N}}.$$

$1/N$ 抵消了选中一个模式时得到的正交和 $N$。这也说明标准 DFT 并没有自动输出“原信号的物理峰值”。

4. 完整手算 4 点 DFT

输入 $x=[1,0,-1,0]$，$N=4$：

$$X[k]=1-e^{-j\pi k}.$$

$k$	$e^{-j\pi k}$	$X[k]$
0	1	0
1	-1	2
2	1	0
3	-1	2

所以 $X=[0,2,0,2]$。逆变换：

$$\begin{aligned} x[n] &=\frac14(2e^{j\pi n/2}+2e^{j3\pi n/2})\\ &=\frac12(e^{j\pi n/2}+e^{-j\pi n/2})\\ &=\cos(\pi n/2), \end{aligned}$$

在整数 $n=0,1,2,3$ 恢复 $[1,0,-1,0]$。$k=3$ 对应负频率 $-\pi/2$，因此两边构成一对共轭正弦分量。

5. FFT 如何减少运算

直接 DFT 有 $N$ 个频点，每点求 $N$ 项，主要计算量为 $O(N^2)$。偶数 $N$ 时，把采样下标分为偶数和奇数：

$$\begin{aligned} X[k] &=\sum_{m=0}^{N/2-1}x[2m]W_N^{2mk} +\sum_{m=0}^{N/2-1}x[2m+1]W_N^{(2m+1)k}\\ &=E[k]+W_N^kO[k]. \end{aligned}$$

因为 $W_N^2=W_{N/2}$，$E$ 与 $O$ 分别是偶数、奇数子序列的 $N/2$ 点 DFT。再利用：

$$W_N^{k+N/2}=-W_N^k,$$

得到另一半：

$$\boxed{ X[k]=E[k]+W_N^kO[k],\quad X[k+N/2]=E[k]-W_N^kO[k] }.$$

一对结果只需同一个中间乘积，称为蝶形。$N=2^p$ 时反复拆到单点：

$$C(N)=2C(N/2)+O(N)=O(N\log_2N).$$

对上例，偶数子序列 $[1,-1]$ 的 2 点 DFT 为 $[0,2]$；奇数子序列 $[0,0]$ 为 $[0,0]$，合并得到 $[0,2,0,2]$。

这推导的是基 2 FFT；其他长度有其他分解与算法，FFT 长度并非必须是 2 的幂。

6. 频率轴、实际点数与补零

采样率 $f_s$，不补零时：

$$f_k=\frac{k f_s}{N},\qquad \Delta f=\frac{f_s}{N}.$$

对于负频率部分，$k>N/2$ 可写成 $(k-N)f_s/N$。偶数长度 $k=N/2$ 位于 Nyquist 点，正负表示等价。

实际数据有 $N$ 点，补零到 $M\ge N$ 点：

$$X_M[k]=\sum_{n=0}^{N-1}x[n]e^{-j2\pi kn/M},$$ $$\Delta f_{\mathrm{grid}}=\frac{f_s}{M},\qquad T_{\mathrm{record}}=\frac N{f_s}.$$

最后一个样本时刻为 $(N-1)/f_s$，频谱记录窗通常按覆盖 $N$ 个采样间隔计算 $N/f_s$。这两个时间描述用途不同。

例如 $f_s=1000$ Hz、$N=1000$、$M=4000$：实际观察仍是约 1 s，频点从 1 Hz 加密到 0.25 Hz。补零读取的是同一段有限数据 DTFT 上更多位置；不会增加新的时间信息，也不自动提高分开相邻正弦的能力。

7. 泄漏从哪里来：直接展开几何和

复正弦 $x[n]=Ae^{j\Omega_0n}$，截取 $N$ 点。第 $k$ 个 DFT 值：

$$X[k]=A\sum_{n=0}^{N-1}e^{j(\Omega_0-2\pi k/N)n}.$$

令 $\theta=\Omega_0-2\pi k/N$，用几何和：

$$X[k]=A\frac{1-e^{jN\theta}}{1-e^{j\theta}}.$$

用 $1-e^{j\beta}=-2je^{j\beta/2}\sin(\beta/2)$：

$$\boxed{ X[k]=Ae^{j(N-1)\theta/2} \frac{\sin(N\theta/2)}{\sin(\theta/2)} }.$$

当 $\theta=0$，取极限得到 $AN$；在恰好整周期、对准频点的情况下，其他不同整数频点处为 0。若频率没有对齐网格，很多频点非零，就表现为泄漏。

对真实余弦，要把正负两支复指数都加上。观察记录窗不等于物理信号永远从记录起点重复；DFT 的周期延拓是表示上的约定。

8. 加窗、相干增益与窗能量

令实际数据点数为 $N$，窗 $w[n]$ 同样有 $N$ 点：

$$v[n]=w[n]x[n].$$

定义：

$$S_1=\sum_{n=0}^{N-1}w[n],\qquad S_2=\sum_{n=0}^{N-1}|w[n]|^2,\qquad CG=\frac{S_1}{N}.$$

$S_1$ 用于相干正弦的幅值校正，$S_2$ 用于功率密度归一化，它们不是同一个量。矩形窗 $w=1$ 时 $S_1=S_2=N$。

例如周期型 Hann 窗，$N\ge4$：

$$w[n]=\frac12-\frac12\cos(2\pi n/N).$$

一个整周期余弦和为 0，因此 $S_1=N/2$。平方：

$$w[n]^2=\frac14-\frac12\cos(2\pi n/N) +\frac14\cos^2(2\pi n/N).$$

余弦平方和为 $N/2$，所以 $S_2=3N/8$。于是 $CG=1/2$。对称型 Hann 用 $N-1$ 作分母，有限长度的这些值会不同；应按实际窗计算，不盲套常数。

窗降低旁瓣通常伴随更宽主瓣。按任务选择窗：分开近邻分量、读取孤立峰幅值、估计宽带噪声，关注点不同。

9. 单边峰值幅度为什么乘 2

矩形窗、整周期余弦、$0<k_0<N/2$：

$$x[n]=A\cos(2\pi k_0n/N+\phi) =\frac A2e^{j\phi}e^{j2\pi k_0n/N} +\frac A2e^{-j\phi}e^{-j2\pi k_0n/N}.$$

由正交性：

$$X[k_0]=\frac{AN}2e^{j\phi},\qquad X[N-k_0]=\frac{AN}2e^{-j\phi}.$$

只保留正频率并表示真实余弦峰值：

$$\boxed{\widehat A[k_0]=\frac{2|X[k_0]|}{N}}.$$

$\arg X[k_0]=\phi$ 是以记录起点和余弦为基准的相位。正弦的相位还需换成相同余弦约定。

加窗时：

$$V[k_0]=\frac A2e^{j\phi}S_1+ \frac A2e^{-j\phi}\sum_nw[n]e^{-j4\pi k_0n/N}.$$

第二项是负频率镜像通过窗产生的贡献。如果该项为零，$2|V[k_0]|/|S_1|$ 精确给出峰值；若镜像足够小，则可近似使用。靠近 DC、Nyquist、邻近峰或偏离频点时，要考虑镜像、泄漏和窗的栅栏损失，不能把简单公式视为任意信号的精确幅值。

DC 与偶数长度的 Nyquist 点各只有一个独立频点，不乘 2。对实信号、$M$ 点 FFT 的单边范围 $0\le k\le\lfloor M/2\rfloor$：

• $k=0$ 不倍增。
• $M$ 为偶数时，$k=M/2$ 不倍增。
• 其余正频率倍增；奇数长度的最高正频率仍要倍增。
• 复数输入没有同样的共轭对称，不能一般地丢弃负频率。

补零后分母仍由 $S_1$ 决定。矩形窗就是实际 $N$，不是补零后的 $M$。

10. 幅值、RMS、功率和 PSD

量	电压输入时的单位	回答
峰值幅度	V	一个正弦峰值是多少
RMS 幅度	V	与同等均方对应的幅度
均方/功率型量	V²	电压平方的平均
PSD	V²/Hz	每单位频率的均方密度
ASD，即 $\sqrt{\mathrm{PSD}}$	V/$\sqrt{\mathrm{Hz}}$	密度的平方根

零均值正弦 $x=A\cos$ 的均方 $A^2/2$，RMS 为 $|A|/\sqrt2$。PSD 的积分得到均方，再开方得到 RMS；单个 PSD 峰高会随窗、实际记录长度和估计方法改变，不能直接当成峰值电压。只增加补零长度时，同一物理频率处的密度值不变；更密的网格可能采到更接近连续谱峰顶的位置，因此网格上的最大值可以变化。这与改变 PSD 归一化是两回事。

电阻上的瓦特还需 $P_{\mathrm W}=x_{\mathrm{RMS}}^2/R$。将 V² 随意标成 W 会导致物理单位错误。

11. 从 Parseval 推导 PSD 归一化

对 $M$ 点 FFT，实际加窗数据 $v[n]$ 在 $n\ge N$ 处补零。由 DFT 正交性：

$$\begin{aligned} \sum_{k=0}^{M-1}|V_M[k]|^2 &=\sum_{n,m}v[n]v[m]^* \sum_k e^{-j2\pi k(n-m)/M}\\ &=M\sum_{n=0}^{N-1}|v[n]|^2. \end{aligned}$$

定义双边周期图密度：

$$\boxed{P_2[k]=\frac{|V_M[k]|^2}{f_s S_2}}, \qquad \Delta f=\frac{f_s}{M}.$$

对全部频点积分近似：

$$\begin{aligned} \sum_{k=0}^{M-1}P_2[k]\Delta f &=\frac1{MS_2}\sum_k|V_M[k]|^2\\ &=\boxed{\frac{\sum_{n=0}^{N-1}|w[n]x[n]|^2}{S_2}}. \end{aligned}$$

矩形窗时右侧正好是原记录均方。一般窗得到窗能量加权的均方；对适当平稳过程具有相应估计解释，并非每段任意确定信号都与未加窗均方精确相等。

这个归一化与 SciPy 的频谱估计说明 及 periodogram 文档 一致。本章不默认去均值；实际库调用可能默认去均值，比较 DC 和总均方前应核对设置。

12. 单边 PSD 为何也是乘 2，但不是幅度平方乘 4

对实数据，负频率与对应正频率的平方模相同。构造单边 PSD 时，把两边的功率相加：

$$P_1[k]= \begin{cases} P_2[k],&k=0\\ P_2[k],&M\ \text{偶数且 }k=M/2\\ 2P_2[k],&\text{其余单边频点}. \end{cases}$$

因此单边 PSD 积分仍等于双边积分。峰值幅度乘 2 与功率合并乘 2 是不同操作；不能先构造单边峰值，再平方直接冒充 PSD，否则还缺时间/窗/采样率归一化。

小例子把单位对上

$x=[1,0,-1,0]$ V、$f_s=4$ Hz、矩形窗，$N=M=4$、$\Delta f=1$ Hz：

$$X=[0,2,0,2].$$

正频率 1 Hz 的单边峰值：

$$A=\frac{2\times2}{4}=1\ \mathrm V.$$

双边密度在 $k=1,3$ 都是：

$$P_2=\frac{2^2}{4\times4}=0.25\ \mathrm{V^2/Hz}.$$

合成单边：

$$P_1(1\,\mathrm{Hz})=0.5\ \mathrm{V^2/Hz}.$$

积分 $0.5\times1=0.5$ V²，RMS 为 $\sqrt{0.5}=0.707107$ V，与时域 $\sqrt{(1+0+1+0)/4}$ 一致。

补零到 $M=16$：网格变成 0.25 Hz，原 1 Hz 点的变换值仍是 2，幅值分母仍是 $N=4$，PSD 分母仍是 $f_sS_2=16$。新增频点让窗频谱曲线更密，但全部单边 PSD 乘 0.25 Hz 后的总和仍为 0.5 V²。

13. 等效噪声带宽与平均

对于实窗且 $S_1\ne0$：

$$\boxed{\mathrm{ENBW}=f_s\frac{S_2}{|S_1|^2}}.$$

它把功率谱型归一化 $|V|^2/|S_1|^2$ 与 PSD 归一化 $|V|^2/(f_sS_2)$ 联系起来。矩形窗 ENBW 为 $f_s/N$；上述周期 Hann 为：

$$f_s\frac{3N/8}{N^2/4}=1.5f_s/N.$$

补零不改变 $S_1,S_2$，因此不改变这个噪声带宽。单次周期图的随机起伏通常很大；Welch 方法分段加窗并平均，减少方差，同时短段会扩大频谱主瓣。段长、重叠和窗应随结果记录。

14. FFT 卷积为何需要补到足够长度

长度 $L_x$ 与 $L_h$ 的线性卷积长度是 $L_x+L_h-1$。DFT 乘法给循环卷积：

$$y_M[n]=\sum_{m=0}^{M-1}x[m]h[(n-m)\bmod M].$$

若 $M<L_x+L_h-1$，尾部会绕回开头，造成时域混叠。两段都补零到：

$$M\ge L_x+L_h-1$$

后，FFT 相乘再逆变换才能得到完整线性卷积。例如 $[1,2,3]*[0.5,0.5]$ 长度为 4；3 点循环卷积会把最后的 1.5 叠回第一点。

15. 实际分析顺序

1. 记录实际 $f_s,N$、信号单位、是否去均值或去趋势。
2. 根据目标选窗，计算实际 $S_1,S_2$。
3. 按需要补零到 $M$，用 $f_s/M$ 建频率轴。
4. 选择峰值/RMS/PSD 的一种表示，按其规则归一化。
5. 实输入构造单边时处理 DC、Nyquist 和长度奇偶。
6. 用已知正弦核幅值，用 Parseval 核总均方。
7. 说明记录长度、窗、平均方式；不要只报告“用了 FFT”。

16. 本章练习与答案

1. $[1,1,1,1]$ 的 DFT：$[4,0,0,0]$；DC 为 $4/4=1$，不乘 2。
2. $[1,-1,1,-1]$ 的 DFT：$[0,0,4,0]$；Nyquist 点幅值为 1，不乘 2。
3. 实际 1000 点、FFT 4000 点，矩形窗幅值分母？1000；频点间隔分母？4000。
4. 单边 PSD 积分为 4 V²，含 DC 的总 RMS？2 V；交流 RMS 需先处理 DC。
5. 5 点实 FFT 的最高正频率点是否倍增？是，奇数长度没有独立 Nyquist 点。