驱动 AY-3-8910（二）：简单频谱分析以及滤波电路设计

前言

在上一篇文章[[驱动 AY-3-8910（一）]]中，已经成功驱动了 AY-3-8910，并让它能够播放指定音量和时长的音调。在这一篇文章中，将介绍频谱分析，以及简单滤波电路设计等内容。

频谱分析

在上一篇文章中，芯片输出信号是直接驱动耳机的，这种未经优化的信号听起来非常刺耳，尤其在高音区。相比之下， 同样是高音，为什么像施坦威这类名琴发出的声音，反而显得清脆悦耳、毫无噪感呢？

要回答这个问题，就得先了解声音的基本特征。钢琴弹奏一个 C4 的音，和小提琴拉奏一个 C4 的音，这两者的音高和音量完全一致，就算蒙住眼睛，也能轻易分辨出是哪个乐器发出的。因为它们的音色 (Timbre) 有显著不同。音色的本质是什么？为什么不同的乐器会有不同的音色？

以钢琴为例，按下一个键时，会联动一个小木槌敲击琴弦振动发声。如果去看钢琴声音的波形图，放大之后就会发现并非简单的正弦波，而是一种更加复杂的振动模式。因为除了琴弦振动本身带来的正弦波以外，还会有诸如琴弦材质、小木槌的摩擦、腔体的共振乃至空气湿度等等因素影响着最终波形的塑造。

在这个过程中，琴弦本身的振动就是基频 (Fundamental Fre­quency)。基频决定了我们听到的音高；除了基频以外的其他所有频率分量就构成了泛音 (Overtone)。在泛音中，基频的整数倍 $2f, 3f, 4f$ 等为谐音 (Harmonics)，这是乐器音色的主要决定部分。除此之外还有一些非谐波的成分，例如小木槌的摩擦、腔体共振等等，这些共同构成了钢琴丰富且复杂的音色。而正是因为不同乐器的发声结构千奇百异，才会有独属于它们自己的音色。

傅里叶变换

知道了构成音色的到底是什么以后，该如何对其进行分析？法国数学家傅里叶 (Joseph Fourier) 揭示了一个深刻的真理：任何复杂的信号波形，本质上都可以分解为无数个不同频率正弦波的线性叠加。这一数学工具便是名垂青史的傅里叶变换 (Fourier Transform, FT)。以钢琴弹奏出的某个音符为例，宏观上看，它表现为一种极度复杂的振动模式；这种复杂的波形通过 FT 后会被分解为一系列正弦分量，背后的频率构成便一目了然。这一原理不仅是声学的基础，更是现代通信、图像处理及信号分析等领域的灵魂。

录音时的电平图像、示波器上显示的波形等，记录的都是物理量随时间流逝的状态。横坐标是时间，纵坐标是物理量的值。这种以时间为参考来衡量动态世界的方式被称为时域 (Time Domain) 分析。而通过傅立叶变换，得以从另一个完全不同的视角去审视信号——这就是频域 (Frequency Domain) 分析。前者好比整首曲子的录音，后者就是这首曲子的乐谱。

接下来的讲解来自于 3Blue1Brown 和 Veritasium 的视频，我非常推荐去观看原视频。此处仅做总结，加上一部分自己的感悟。

从数学上讲，连续信号的傅里叶变换通过积分将时域函数 $f(t)$ 映射为频域函数 $F(\omega)$：

$$F(\omega) = \int_{-\infty}^{\infty}f(t)e^{-j\omega t}\mathrm{d}t$$

这个公式乍一看很吓人，不妨拆开来看看。公式中的 $j$ 就是复数单位 $i$. 在复平面中，横坐标代表实部，纵坐标代表虚部。这就意味着当一个矢量乘以 $i$ 后方向会和原矢量垂直。

假设某一函数反映了平面内某一点的位置随时间变化的关系：

$$\frac{ \text{d} }{ \text{d}t }e^{it}=i\cdot e^{it}$$

把等式左侧看作速度矢量，右侧看作位置矢量。由于右侧比左侧多了一个旋转因子 $i$，速度矢量的方向将会永远垂直于位置矢量。随着时间 $t$ 的增长，将在复平面上画出一个圆。这便引出了数学中最伟大的公式——欧拉公式 (Euler's Formula)：

$$e^{it} = \cos(t)+i\sin(t)$$

如果要改变绘制这个圆的速度，需引入频率参数 $\omega$ 改变速度矢量的虚部 $\omega i\cdot t$。如果再为 $e$ 的幂加上实部变成 $(\sigma+\omega i)t$，此时速度矢量在水平方向上也有了分量，随时间变化就能够绘制出螺旋线。

将旋转基底乘以原始信号函数 $f(t)$，在视觉上这就产生了一种类似于将原信号「缠绕」在单位圆上的效果。在这个缠绕图像中，质心的位移正是频率提取的关键。从物理上看，质心的位置取决于图像分布的均衡程度。计算质心的方法是取若干个样本点，并计算它们位置的平均值。而当样本点的个数趋于无穷时，就变成了定积分。当缠绕频率与信号内在频率相同时，质心会猛然偏离原点，偏移距离和方向便反映了该频率的强度与相位。

为了方便解释，故采用了质心这一思想。而傅里叶变换就不需要对图形的每个采样点取平均值。换个说法，最大位移就指向了原函数的频率。

傅里叶变换要求信号必须是收敛的，即必须衰减得足够快，在 $t\to\infty$ 内能够收敛到 0 或一个常数上，以及能量是有限的。换句话说，信号是绝对可积的：

$$\int_{-\infty}^{\infty}|f(t)|dt < \infty$$

以及平方可积的：

$$\int_{-\infty}^{\infty}|f(t)|^2dt < \infty$$

而如果信号是永恒不变的，在尝试计算前文提到的质心的时候，积分的结果就会发散，从而无法反映出正确的结果。直流和阶跃是两种很重要的信号，然而它们并不满足绝对可积的条件。这时候就要引入一个衰减因子，也就是前文中提到的在 $e$ 的幂中的实部 $\sigma$，强行让信号收敛。而这就是拉普拉斯变换 (Laplace Transform)：

$$F(s)=\int_{0}^{\infty}f(t)e^{-st}\text{d}t$$

其中，拉普拉斯算子 $s = \sigma + j\omega$. FT 可以近似地看作 LT 在特定条件下的一种特殊情况。

在工程中通常不会对时间跨度这么长的信号进行分析，而是从中截取一段有限的时间区间 $[0, T]$ 进行分析，就像只对屏幕上窗口内展示的信号进行分析，原信号 $f(t)$ 乘以一个窗函数 (Window Function)。在这段区间以外的信号都会强制变成 0.

然而，这种截断并非没有代价。如果截取的区间未能恰好包含信号频率的完整周期时，窗口的首尾两端就会发生剧烈的跳变。这种时域上的不连续性就会造成很多本不应该出现的频率分量，导致频谱能量从主频向四周漫延。这便是所谓的频谱泄露 (Spectral Leakage)。频率主峰附近会隆起许多本不属于信号的旁瓣 (Sidelobes)。

为了抑制这种副作用，信号应该「淡入淡出」。通过使用加权窗，窗口中心区域的信号获得更高的权重，而让两侧边缘的信号平滑地衰减至零。这种柔和的过渡能够有效平复截断处的跳变，从而压制旁瓣能量。常见的窗函数如下，此处暂时不做过多解释：

读者可能已经注意到了，频率相同时的图形颇有意思，此处粗略介绍一下。假设信号是一个简谐波加上直流偏量：$g(t)=B+A\sin(2\pi f_0t)$. 缠绕的过程本质上是将时间 $t$ 映射为极坐标中的角度 $\theta$. 如果以频率 $f$ 进行旋转，那么在 $t$ 时刻，旋转过的总角度（弧度）为：

$$\theta = 2\pi ft$$

可得时间 $t$ 与角度 $\theta$ 的关系：

$$t=\frac{\theta}{2\pi f}$$

代入 $g(t)$ 作为极径 $r$：

$$\begin{align} r(\theta) &= g\left( \frac{\theta}{2\pi f} \right)=B+A\sin\left( 2\pi f_0\cdot \frac{\theta}{2\pi f} \right) \\ &= B+Asin(\frac{f}{f_0}​​\theta) \end{align}$$

而当 $f=f_0$ 时就会变成：

$$r(θ)=B+Asin(θ)$$

这正是帕斯卡蜗线 (Limaçon of Pascal) 的标准方程。而当 $B=A$ 时就是大名鼎鼎的心脏线 (Cardioid)。

回到正题。现实生活中各种信号变化并非无限精细，同时仪器输出的波形看似光滑，也会存在不可分辨的最小单位。所以最终得到的是离散、有限的信号。此时连续傅里叶变换不再适用，取而代之的是离散傅里叶变换 (Discrete Fourier Transform, DFT)。对于长度为 $N$ 的序列 $x[n]$，对其 DFT 的公式为：

$$X[k] = \sum_{n=0}^{N-1}x[n]\cdot e^{-j\frac{2\pi}{N}nk}, k=0, 1, 2, \dots, N-1$$

离散的视角下，圆周不能平滑地旋转，而是像那种一次走一格的时钟。在 $e$ 的幂中，如果有 $N$ 个采样点，那么就需要将圆弧分成 $N$ 份，也就是 $\tfrac{2\pi}{N}$. 而 $n$ 就表示当前是第几个采样点。$k$ 就是缠绕频率。

DFT 可以看作一个矩阵乘法 $y=Wx$，其中 $W$ 是 $n\times n$ 的范德蒙德矩阵：

$$\begin{bmatrix} X[0] \\ X[1] \\ \vdots \\ X[n-1] \end{bmatrix} = \begin{bmatrix} \omega^0 & \omega^0 & \dots & \omega^0 \\ \omega^0 & \omega^1 & \dots & \omega^{n-1} \\ \vdots & \vdots & \ddots & \vdots \\ \omega^0 & \omega^{n-1} & \dots & \omega^{(n-1)(n-1)} \end{bmatrix} \begin{bmatrix} x[0] \\ x[1] \\ \vdots \\ x[n-1] \end{bmatrix}$$

为了得到单个点 $X[k]$，需要从 0 到 N-1 进行求和，其中涉及到 N 次复数乘法与 N-1 次复数加法，复杂度就是 $O(n)$；而整个序列中含有 N 个点，需要将上述过程重复 N 次，总复杂度就是 $O(N^2)$. 对于高采样率的信号计算量会变得非常巨大，无法实时分析。

快速傅里叶变换 (Fast Fourier Transform, FFT) 则优化了算法的复杂度。这个算法曾经被数学家吉尔伯特·斯特朗描述为「我们一生中最重要的数值算法」。库利-图基算法是最常见的 FFT 算法。除此之外还有质因子算法 (Prime-Factor Algorithm) 等等，此处按下不表。

库利-图基算法靠两个手段优化。第一个是旋转因子的周期性与对称性。由于旋转因子 $W^{nk}_N = e^{-j\frac{2\pi}{N}nk}$ 本质上就是圆上的点，因此具备两个特性：对称性 $W^{k+N/2}_{N} = -W^{k}_{N}$ 和周期性 $W^{k+N}_{N} = W^k_N$. 这意味着很多复杂的乘法运算结果其实是一样的，或者只需变个符号。

而第二个是递归拆分。将一个长度为 $N$ 的大问题，拆分成两个长度为 $N/2$ 的小问题。对于长度为 $N$ 的序列 $x[n]$，分成偶数索引 $x_{2m} = x_0, x_2, \dots, x_{N-2}$ 和奇数索引 $x_{2m+1} = x_1, x_3, \dots, x_{N-1}$ 两组。分别计算这两组的 DFT。 利用旋转因子的特性，将这两个小结果组合成最终的结果。这个过程会一直递归下去，直到拆成只有 2 个点的小单元。此时算法的复杂度就降到了 $O(N \log_2 N)$.

得益于 FFT 的高效，它不仅统治了频率分析领域，更成为了数字信号处理的核心底座。例如，计算机中极为耗时的卷积 (Convolution) 运算，在频域下只需进行简单的乘法即可完成。这种运算在图像、音频处理，以及深度学习（例如卷积神经网络）等领域都有大量运用，以后有机会深入研究。

滤波器

前面的频谱分析结果揭示了方波在高频处有频率分量的分布，而这也是导致听感刺耳的「罪魁祸首」。因此，需要将这些高频的成分过滤掉，实现这个功能的电路就被称为滤波器 (Filter)。滤波器有很多种，如果根据频率响应来划分，通低频、阻高频的就是低通滤波器 (Low-pass Filter, LPF)，反之就是高通滤波器 (High-pass Filter, HPF)；还有针对特定频段的带通 (Band-pass Filter, BPF) 和带阻 (Band-stop Filter, BSF) 滤波器。

除此之外，有些滤波器结构比较简单，可能只包含电阻、电容或电感这几种器件，仅起滤波作用，这类就被称为无源滤波器；而有的含有运算放大器、晶体管等元件，除了滤波还可以对信号提供增益，这类就被称为有源滤波器。在本篇文章中，将着重介绍最简单的 RC 无源低通滤波器。

这个电路的结构非常简单。如图所示，将一个电阻 R1 与信号输出串联，电容器 C1 与负载 RL 并联，就可以产生 RC 低通响应。电容器本质上就是由两块导体板与中间的绝缘介质（电介质）构成的。在直流电路中，当电容器接入电路中时，电极的电压会因为电子的聚集而不断上升，直到等于电源电压，电流停止流动，所以此时表现可近似看作「断路」。换句话说，对电流的阻抗 (Impedance)，或容抗 (Capacitive reactance) 很大。

而在交流电路中，由于电流方向不断反转，电压方向正向时，电容板充电；当电压方向反向时，电容板放电，并向相反方向充电。电流并没有实际流过绝缘介质，但电流在电路中不断地来回流动，就像在电路中「通过」了一样；而频率越高，电容器充放电的速度也越快，内部介质来不及建立稳态电压，电流在外部电路中来回流动的速率也越高。从电路外部来看，这等效于电容器的阻抗越小。总结来说，电容的一个重要性质就是「通高频，阻低频」。容抗的计算公式如下，也可以看出其大小与频率成反比：

$$W X_C = \frac{1}{2\pi fC}=\frac{1}{\omega C}$$

当输入低频信号时，容抗相对于电阻的阻抗更高，因此电容上的分压降低；当输入高频信号时，容抗相对于电阻的较低，这意味着电阻上的分压降低，这样就较少的电压传输到负载。因此，低频通高频阻。而反过来，如果想要实现高通滤波，R1 和 C1 的位置就需要对调一下。这样，低频分量就会更多地分压在 C1 上，而高频分量更多地分压在 R1 上。

额外补充一点，除了电容器以外，电感器 (Inductor) 也有着与之相似的电气特性，只是两者恰好相反，即「通低频阻高频」。从物理本质上看，电感器通常由导线绕制成的线圈构成。当电流流经线圈时，会在其周围激发出磁场。如果电流 $i$ 的大小或方向发生变化，穿过线圈的磁通量 $\Phi$ 也会随之改变。根据法拉第电磁感应定律 (Faraday's Law of Induction)，当线圈不动而磁场变化时，穿过回路的磁通量也发生变化，由此在回路中激发的感生电动势 $\mathscr{E}$ 的大小与穿过线圈的磁通量变化率成正比。其数学表达式为：

$$\mathscr{E} = -n\frac{\mathrm{d}\Phi}{\mathrm{d}t}$$

其中，$n$ 是线圈的匝数。公式中的负号就是楞次定理 (Lenz's Law) 的体现，它表明感生电动势的方向始终倾向于阻碍原磁通量的变化。这种性质在电路中表现为对交流信号的阻碍作用，称为感抗 (Inductive Reactance)。感抗的计算公式如下：

$$X_L = 2\pi fL=\omega L$$

可以看出，感抗的大小与信号频率 $f$ 成正比。在低频或直流（$f=0$）情况下，感抗极小，电感可近似看作导线；而在高频情况下，感抗显著增大，从而实现对高频噪声的有效滤除。

回到正题。「高频」和「低频」是个相当宽泛的概念。要落地到电路设计，必须计算出具体数值。滤波器不会引起显著衰减的频率范围称为通带，反之则称为阻带。模拟滤波器（如 RC 低通滤波器）在滤波并不会像一堵墙一样瞬间把信号切断，而是像一个陡坡，让信号强度从通带逐渐衰减到阻带。这意味着我们很难在曲线上找到一个「绝对」的阻塞起点。在工程实践中，通常将信号强度衰减至输入端 -3 dB 处的频率定义为截止频率 (Cutoff Frequency)。这里使用分贝 (Decibel, dB) 来定量描述信号功率的损耗，这种对数单位在衡量系统增益或衰减时比线性数值更加直观。为什么偏偏是 -3 dB 呢？由这个单位的计算方法就能得出：

$$L_{\mathrm{dB}} = 10 \log_{10} {\left( \frac{P_1}{P_0} \right)}$$

$P_0$ 是输入功率，$P_1$ 是输出功率。当 $\tfrac{P_1}{P_0} = \tfrac{1}{2}$ 时，$L_{\mathrm{dB}} \approx -3.01 \text{ dB}$，正好对应信号功率衰减到最大功率的一半。在该点上，电阻的阻值等于电容的容抗，即：

$$R = X_C$$

代入容抗的计算公式，即可得到截止频率：

$$f = \frac{1}{2\pi RC}$$

信号采样与分析

为了对芯片输出的音频信号进行频域分析，需要准确采样信号。这里我用一块 USB 声卡配合 Adobe Audition 和示波器进行录制与分析。

操作时有个关键的技术细节，信号必须接在声卡的 LINE IN（线路输入）接口上，而不是 MIC IN（麦克风输入），因为两者的电平定义完全不同。

与带有预放大功能的 MIC IN 不同，LINE IN 专为高电平信号设计。大多数的 PSG 芯片输出信号都会直接用于驱动小型的模拟电路或扬声器，它们通常在 0 V 到 5 V 之间做单极性跳变，峰峰值 $V_{p-p}$ 能达到 2 V 甚至 4 V。然而，消费级声卡预期的输入电平应该是 -10 dBV。这个单位是以 1 V 为基准参考的电压值。根据换算公式：

$$\text{dBV} = 20\cdot\log_{10}\left(\frac{V}{1.000}\right)$$

可以计算出对应的有效值电压约为 $0.316\text{ V}_\text{RMS}$.

$\text{ V}_\text{RMS}$ 指的是均方根 (Root Mean Square) 值，是衡量随时间变化的电压有效性的核心指标。比如正弦波为的有效值约等于峰值电压的 0.707 倍（即 $\tfrac{1}{ \sqrt{2} }$）；而理想方波信号的有效值就等于它的峰值。由此可见，芯片输出的原始信号强度远超过声卡的额定承受范围，直接输入会导致严重削波失真 (Clipping)，也就是在电平表上看到的「爆音」。

为了安全且完整地采集信号，需要设计一个简单的调理电路。芯片输出通常携带 2.5 V 左右的直流偏置，这会偏置声卡的输入级。首先串联一个电容器，它在电路中起到了高通滤波器的作用，只允许变化的音频成分通过，阻断直流成分。

随后再引入大电阻分压网络来进一步衰减信号强度，将 5 V 的 TTL 电平降压至 0.3 V 左右。这里就涉及到阻抗匹配 (Impedance Matching) 的概念。信号传输可以简化为一个串联分压模型：

$$V_{in} = V_{src}\cdot\frac{ Z_{in} }{ Z_{out}+Z_{in} }$$

其中，$V_{src}$ 和 $Z_{out}$ 为信号内部的理想电压和输出端的阻抗（源阻抗），而 $Z_{in}$ 和 $V_{in}$ 为输入端的阻抗（负载阻抗）和声卡实际接收到的电压。为了让声卡尽可能完整的接收到芯片产生的电压信号，即让 $V_{in}$ 尽可能接近 $V_{src}$，那么 $Z_{in}$ 应该尽可能大。当 $Z_{in}\to+\infty$ 时，$V_{in}=V_{src}$. 通俗点说，输出阻抗越低，信号源的「推力」就越强，信号也就越不容易在传输过程中损耗。

为了观察滤波器的实际效果，先录制一段未处理的音符 C4（注意音量大小）：

note.PlayNoteEvent(500, NOTE_C4, 15);

随后，在信号输出端与地之间并联了一个 10nF 的陶瓷电容。这个元件与前级的限流电阻共同构成了一个一阶低通滤波器。在示波器上可以明显观察到，原本陡峭的方波边缘变得圆润了——这意味着高频谐波被削弱，信号在频域上的能量分布发生了显著变化。经过滤波的音频听上去要更「闷」一点。

由于输出端的内阻未知，可以通过逆向工程来获取系统的关键参数。利用示波器捕捉信号的阶跃响应。实测显示，信号的上升沿时间（10%-90%）$t_r$ 为 40 μs 左右。由于示波器的输入阻抗很大（1 MΩ），远大于芯片的输出内阻，其并联效应对波形的影响可以忽略不计。对于一阶 RC 滤波器，上升时间 $t_r$ 与时间常数 $\tau=RC$ 满足如下关系：

$$t_r = 2.2\cdot \tau$$

可得芯片的输出内阻约为 1.82 kΩ. 再代入截止频率的计算公式 $f_c=\tfrac{1}{2\pi RC}$ 可以推导出：

$$f_c\approx\frac{0.35}{t_r}$$

得截止频率 $f_c\approx8.75\text{ kHz}$. 这个结果和 Audition 中频谱分析展示的结果非常接近。可以看到，信号在 8 kHz 附近恰好产生了约 -3 dB 的衰减。

受限于奈奎斯特采样定律 (Nyquist Sampling Theorem)，以常见的音频采样率 44.1 kHz 为准，声卡中的抗混叠滤波器会强制切除 22 kHz 以上的所有能量。因此在接近 20 kHz 时开始急速衰减，并在 22 kHz 之后衰减到 -100 dB 以下。

在频域的视角下，方波本质上是由无数正弦波叠加而成的复合信号。对于一个振幅为 $A$，频率为 $f_0$ 的理想方波，其傅里叶级数展开式为：

$$f(t)=\frac{4A}{\pi}\sum_{n=1, 3, 5, \dots}^{\infty}\frac{1}{n}\sin(2\pi nf_0t)$$

这个公式揭示了方波的三个本质：

1. 只包含奇次谐波：方波只包含基波（$f_0$）及其奇数倍频率（$3f_0, 5f_0, 7f_0 \dots$）；
2. 振幅衰减：谐波次数越高，能量越弱。第 $n$ 次谐波的振幅是基波的 $1/n$；
3. 相位对齐：所有的正弦波在零点同时出发，通过叠加抵消，最终在时域上构建出平整的顶部和陡峭的边缘。

在试图用有限项谐波拟合这一理想模型时，在跳变的边缘处总会出现一圈无法消除的振铃。这种现象被称为吉布斯现象 (Gibbs Phenomenon)。即便将谐波项数增加到更多，边缘过冲的幅度依然会锁定在 9% 左右。

现实世界中不存在绝对完美的方波。受限于系统的带宽限制，高频成分在传输中必然被截断，导致波形转角处不可避免地出现圆角或过冲。这也是为什么示波器捕捉到的信号，永远无法呈现出理论模型中那种绝对锐利的直角。

总结

这篇文章主要侧重于理论和电路上的知识。讲的东西已经够多了，从傅里叶变换到滤波器设计，最后再到信号采样与分析，虽然只是稍作介绍，途中也学到了不少新知识。此后将开始尝试解析 Tracker Music 的文件格式，并写出能够播放它们的驱动程序。对于一些有更多轨道的 chiptune，可能会考虑使用性能更强大的单片机，如 ESP32 等。

扩展阅读

1. Fourier transform. Wikipedia. https://en.wikipedia.org/wiki/Fourier_transform
2. Fourier series. Wikipedia. https://en.wikipedia.org/wiki/Fourier_series
3. Fast fourier transform. Wikipedia. https://en.wikipedia.org/wiki/Fast_fourier_transform
4. Window function. Wikipedia. https://en.wikipedia.org/wiki/Window_function
5. Spectral leakage. Wikipedia. https://en.wikipedia.org/wiki/Spectral_leakage