现在是时候把理论付诸实践,制作一些音频滤波器和均衡器了。你知道一个滤波器的系数决定了它的频率响应和其他特性。但你如何找到这些系数呢?有两种常用的方法来计算IIR滤波器的系数。
直接Z面设计
模拟滤波器到数字滤波器的转换
本章使用了以下滤波器的命名规则。
LPF:低通滤波器
HPF:高通滤波器
BPF:带通滤波器
BSF:带阻滤波器
APF:全通滤波器
HSF:高架滤波器
LSF:低搁置滤波器
PEQ: 参量EQ滤波器
GEQ: 图形EQ系统
11.1 直接z平面设计
在这第一类设计技术中,你直接在z平面上操作极点和零点,以创建你想要的响应。你可以利用将系数与极点/零点位置联系起来的简单方程。对于一般的双曲线,我们可以用方程11.1和11.2将系数和极点/零点位置联系起来。
对于分子
对于分母
对于分子或分母,a1或b1系数直接控制着零点或极点的角度。与零点或极点的距离R由a1、a2或b1、b2决定。对于一阶滤波器,系数只控制实轴上的极点和零点的位置。没有共轭对。然而,仔细放置极点和零点仍会影响滤波器的响应。
11.1.1 简单的谐振器
谐振器是一种可以使其具有非常窄的峰值的BPF。基本的谐振器在二阶反馈拓扑结构中实现了一对共轭的极点。图11.1a显示了使用直接形式结构的滤波器的框图。
图11.1: 谐振器的(a)直接形式框图和(b)极点/系数关系。
谐振器的差分方程为:
谐振器的工作原理是对二阶反馈网络形成的共轭极点进行简单操作。b2项控制到极点的距离,使谐振峰更尖锐(当极点接近单位圆时)或更宽(当极点远离单位圆时)。b1和b2项控制极点的角度,控制中心频率,如图11.1b所示。a0项用于缩放滤波器,所以它的峰值输出总是归一增益或0dB。
这种设计的缺点是,频率响应只在一个频率上是对称的,即π/2,此时低频和高频的幅度向量是对称的。在所有其他频率上,响应是不对称的。图11.2显示了谐振器的极点在单位圆的右半部(a)和左半部(b),这导致响应分别向低频(a)和高频(b)倾斜。
图11.2:(a)当谐振器极点在单位圆的右半部分时,响应是不对称的,并向低频一侧倾斜。(b) 当谐振器极点在单位圆的左半部分时,响应是不对称的,并向高频侧倾斜。(频率轴是线性的,以更好地显示不对称性。)
11.1.2 史密斯-安盖尔谐振器
改进后的设计,被称为Smith-Angell谐振器,通过在滤波器中加入两个零--一个在z=1,一个在z=1--来减少不对称问题,以便用零的幅度响应来锁定直流和奈奎斯特。这迫使滤波器在某种程度上变得更加对称(但不是完全对称),并具有使带通在本质上更具选择性的优势。这个设计也是用a0进行增益归一化。和以前一样,先用b2设定半径,然后用b2的系数和所需的极点频率计算b1。滤波器并不是真正的归一化为0.0dB,实际上有一个小于1dB的轻微波动。
Smith-Angell谐振器的差分方程是
图11.3a显示了使用直接形式结构的框图,11.3b显示了极点/零点关系。图11.4a显示了当极点在单位圆的右半部分时谐振器频率响应的形状,而图11.4b将极点移到单位圆的右半部分。
图11.3:Smith-Angell谐振器在不同中心频率设置下的频率响应,(a) fc = 2 kHz,(b) fc = 20 kHz。(对于这两张图,fs = 44.1 kHz;频率轴是线性的,以便与简单的谐振器进行比较)。
图11.4: Smith-Angell谐振器在不同中心频率设置下的频率响应,(a) fc = 2 kHz,(b) fc = 20 kHz。(对于这两张图,fs = 44.1 kHz;频率轴是线性的,以便与简单的谐振器进行比较)。
11.2 模拟滤波器到数字滤波器的转换
一个比较广泛使用的滤波器设计方法是,首先从经典的模拟设计开始,然后将它们转换成数字版本。这样做是有道理的,因为有很多优秀的模拟设计已经为我们完成了。我们只是需要一种方法来使它们在我们的数字世界中工作。虽然模拟滤波器的设计不在本书的范围之内,但这两个设计世界之间有许多相似之处。例如,模拟和数字滤波器的设计都涉及到一个传递函数,你可以通过操作来产生极点和零点。它们也都使用一个变换来从时域到频域。一个根本的区别是,在模拟中,没有奈奎斯特的限制,从-∞到+∞的所有频率都包括在内。模拟滤波器使用电感和电容等反应性元件来完成模拟信号的积分或微分功能。数字滤波器则使用延迟元件来产生相位偏移。表11.1总结了这些异同点。你可以在《用C++设计软件合成器插件》中找到对这两种技术更详细的分析。
在图11.5中,你可以看到右边的s平面,它也是一个复数平面。实轴被称为σ轴,虚轴是jω轴。jω轴是频率轴,它的跨度是-∞到+∞拉德/秒。单位圆映射到jω轴上-Nyquist和+Nyquist之间的部分。为了将模拟滤波器设计转化为数字滤波器设计,我们需要一个映射系统来迫使s面的极点和零点变成z面的极点和零点。一旦我们有了这个系统,就可以计算在数字位置产生这些极点和零点的系数。在模拟世界中,接近jω轴的极点将成为数字世界中接近单位圆的极点。
表11.1: 模拟和数字滤波器设计属性摘要
数字滤波器设计 | 模拟滤波器设计 |
---|---|
使用一个复数域传递函数来联系输入/输出 | 使用一个复数域的传递函数来联系输入/输出 |
延迟元件产生相移和延迟 | 反应式组件进行整合 |
使用zTransform(采样时间到频率) | 采用拉普拉斯变换(连续时间到频率) |
在Z平面上的极点和零点 | S面的极点和零点 |
奈奎斯特限制 | 允许从-oo到+无穷 的所有频率 |
极点 必须在单位圆内才能稳定运行 | 极点必须在s平面的左手部分才能稳定运行 |
图11.5:Z平面和S平面以不同方式显示类似的信息。
问题是,模拟极点和零点可以存在于jω轴的任何地方,甚至在数字奈奎斯特区之外的频率。模拟设计在σ和/或jω轴上的∞和+∞处有一个极点或零点也很常见。在图11.6中,你可以看到这个问题--在½奈奎斯特处靠近jω轴的一对模拟零点将映射到½奈奎斯特角处靠近单位圆的数字零点位置。那么,在s平面上位于奈奎斯特频率之外的一对极点呢,它们应该被映射到哪里?另一种思考这个问题的方式是,在模拟s平面中,整个左手边,包括所有的无限频率,必须映射到图11.7a所示的z平面的单位圆内部。此外,整个右侧平面必须映射到单位圆的外部,如图11.7b所示
因此,我们需要的是一个数学运算来实现从s面到z面的转换。它必须将s平面左侧的无限大面积映射到单位圆在z平面上的有限内部面积,并将s平面右侧的无限大面积映射到无限大面积上
图11.6:模拟极点从s面到z面的映射。
图11.7:要将模拟滤波器转换为数字滤波器,(a)无限大的s平面左侧必须映射到z平面的单位圆内的有限空间,而(b)s平面的右侧平面映射到z平面的单位圆的外部。
在Z平面的单位圆之外。此外,它应该将s平面jω轴上的点映射到z平面上从DC到±Nyquist的单位圆的外缘的相同点。双线性变换可以作为映射装置,将模拟传输函数H(s)转换成数字版本H(z)。有几种不同的方法来推导双线性变换。在这里,我提出一个简单的推导,它依赖于寻找s和z变量之间的关系。你可以在《用C++设计软件合成器插件》中找到一个完全不同的方法,在那里我用一个数字积分器来近似模拟积分的操作。这两种方法都是有效的,尽管第二种方法被认为更加稳健。
如果我们希望将模拟传递函数H(s)转换为采样的等价物H(z),我们需要一种方法来对连续的s面进行采样,以产生一个采样的z面版本。换句话说,我们需要创建一个采样的模拟传递函数Hs(s),其中下标s代表 "采样的"。我们试图找到一个函数g(z),使公式11.7成立。
由于采样间隔是T,那么术语ejωT将对应于一个时间的采样。因此,如果我们通过让z = ejωT来评估H(z),那么我们将得到方程11.8中的采样模拟传递函数Hs(s)。
为了求解,我们注意到,s=jω,z=ejωT,因此。
我们现在有了s和z之间的关系,但是我们需要处理取z的自然对数的问题,这没有一个封闭式的解。如果我们使用ln()的泰勒级数展开,就可以得到Re(z)>0的公式11.10。
这个一般解决方案的第一项近似值就是双线性变换。
我们使用的双线性变换是。
2/T项在后面的数学中会被冲掉,所以我们忽略它。方程11.11通过在无限频率上取值并将其映射到奈奎斯特来执行映射。在-∞或+∞的极点或零点被映射到奈奎斯特的极点或零点。在奈奎斯特和-∞或+∞赫兹之间的其他频率被挤压到奈奎斯特周围的小区域,就在单位圆内,如图11.8所示。
图11.8:双线性变换将左侧平面上奈奎斯特区以外的区域映射到奈奎斯特附近的圆内的区域。
双线性变换通过tan函数将模拟频率非线性地映射到其相应的数字频率,如公式11.12所示。
ωa = 模拟频率
ωd =映射的数字频率
T = 采样周期
在ω的低值时,tan函数是线性的,但随着频率的增加,它变得更加非线性。在低频时,模拟和数字频率是线性映射的。在高频率下,数字频率变得扭曲,不能正确地映射到模拟对应的频率。这意味着,一个给定的具有截止频率fc的模拟设计,在转换后通常会有错误的数字截止频率。解决办法是对模拟滤波器进行预扭曲,使其截止频率在模拟域中处于错误的位置,但在数字域中却处于正确的位置。为了预热模拟滤波器,你只需使用公式11.12应用于模拟滤波器的截止频率。当你把所有的操作结合起来,你就得到了双线性Z变换(BZT)。转换的步骤方法如下。
从具有归一化传递函数H(s)的模拟滤波器开始,其中 "归一化 "意味着截止频率被设定为ω= 1rad/sec;这是描述和指定模拟传递函数的典型方法。
对于LPF和HPF类型,选择数字滤波器的截止频率ωd,用公式11.13得到预扭曲的模拟角频率。
对于BPF和BSF类型,指定ωdL和ωdH(数字滤波器的下角和上角频率),用公式11.14计算两个模拟角频率ωaL和ωaH以及ωo2和W项。
将归一化滤波器的截止频率扩展到新的模拟截止点ωa,方法是在模拟传递函数中用公式11.15中适当的比例系数代替s。
应用双线性变换,用公式11.16代替s。
用代数方法处理数字传递函数H(z),使其变成熟悉的形式,这样你就可以识别系数--这往往是最困难的部分。
11.2.1 示例
将图11.9a中的基本RC模拟LPF转换成数字LPF。采样率为44.1kHz,所需的截止频率为1kHz。
图11.9:(a)一个简单的模拟RC滤波电路。(b) 其数字等效物。
第1步:得到归一化模拟传递函数H(s)。
第2步:计算预扭曲的模拟截断。
第3步:用适当的系数对模拟传递函数H(s)进行去正常化。
第四步:应用BZT。注意最后一步的代数操作,除以z产生z-1项,然后归一化,使分母的第一项为1。
进入标准的H(z)格式。
方程11.20的格式是我们需要的,分子和分母正确组成,以观察系数。从传递函数中可以看出。
a0 = 0.0667
a1 = 0.0667
b1 = -0.8667
差分方程为。
图11.9b所示的直接形式框图显示了一个极点-零的滤波器。这是有道理的:原来的模拟滤波器在s = -1处有一个极点,在无穷远处有一个零点。数字等价物在z = 0.8667处有一个极点,在z = -1处有一个零点(Nyquist),都在实轴上。
11.3 音频Biquad滤波器的设计
在上一个例子中,我们取了一个固定的模拟滤波器,应用BZT将其转换为双线性变换的数字等价物。然而,我们想用同样的概念将具有任意属性的模拟滤波器,如fc、Q或滤波器增益(提升或削减),转换成它们的数字等效物。此外,作为音频工程师,我们还希望这些参数可以连续调节,以便用户可以用旋钮或滑块等GUI元素来控制它们。我们还希望能够轻松地以编程方式调制参数--例如,我们可能想用一个低频振荡器(LFO)来平滑地调制截止频率,使其随时间变化。
好消息是,我们可以应用一些相当直接的三角函数和代数来实现这一目标。此外,有多种方法可以从数学上解决这个问题,所以你可能会发现两套不同的方程,产生相同的频率响应特性。也有多种拓扑结构来创建相同的滤波器(见家庭作业问题4),所以你可能会找到不同的算法,以不同的方式组合biquad结构,但也产生相同的结果。
坏消息是,这个数学问题并不简单:需要做一些工作,包括使用三角函数同义词和一些复杂的代数。幸运的是,DSP工程师在开发信号处理的高效音频算法时,很早就开始研究这个转换问题。其结果是为几乎所有经典的模拟滤波器类型提供了一套强大的滤波器设计方程。如果你想了解更多,可以看看Bristow-Johnson的优秀而简洁的论文(Bristow-Johnson, n.d.)或Motorola的应用说明(见Motorola, 1991),这些都是关于这个问题的最早的工作--你可能想先复习一下三角函数的知识。
以下经典的模拟滤波器被转换为数字形式,并以双通道拓扑结构实现。
LPF
HPF
BPF
BSF
LSF
HSF
恒定和非恒定Q参数滤波器
二阶巴特沃斯LPF和HPF
二阶Linkwitz-Riley LPF和HPF(适用于分频器)
1阶和2阶APF
11.3.1 音频Biquad滤波器的结构
图11.10:音频biquad结构允许干信号和处理后的信号混合。这个版本使用的是直接形式(如虚线框内所示),但任何一个biquad结构都可以代替它。
对于音频专用的滤波器,我们将使用修改过的biquad结构,使我们能够以任何比例混合干的(未受影响的)音频信号和湿的(处理过的)部分。对于许多滤波器的设计来说,混合将始终是相同的:0%的干信号和100%的处理信号。我们将在最终的拓扑结构中加入两个新的系数d0(干)和c0(处理)。你可以选择第10章中的任何一个biquad结构来实现信号处理部分。图11.10显示的是音频双阙结构。你可以看到干、湿信号路径是如何与添加的系数混合的。差分方程为:。
11.3.2 经典滤波器
下面这组滤波器设计实现了基本的模拟滤波器传递函数,这些函数从模拟信号处理理论的早期就已经存在了。有时它们被称为经典滤波器。在所有这些滤波器算法中,输出是100%的处理和0%的干燥,因此在每种情况下。
由于模拟滤波器转换中使用的三角函数,以及BZT需要tan函数的事实,计算通常需要某种三角函数:sin、cos和/或tan。在这些函数中,tan函数有一个问题:它在π/2处没有定义。当设计方程需要调用tan函数时,我们应该总是检查参数以避免未定义条件。我们可以把参数的值夹在非常接近的地方,如0.95π/2。你将在第11.6节中看到在AudioFilter对象中检查和夹紧这个值的代码。这些滤波器在从直流到略低于奈奎斯特的范围内都很稳定。一些设计计算需要一个tan(πfc/fs)的项,当fc=Nyquist时,它将是未定义的。
对于我们的插件控制,我们将始终使用完整的、10倍频程范围的截止或中心频率,在20.0Hz到20480Hz的10倍频程范围内,采样率为44.1kHz;在这个采样率下,20480Hz的上限大约是奈奎斯特频率的93%。如果你的滤波器在这个范围内不稳定,请仔细检查系数的计算,小心tan函数的问题。
11.3.2.1 一阶LPF和HPF
图11.11:(a)一阶LPF和(b)fc=100 Hz、1000 Hz和10 kHz的HPF。
图11.11a和b显示了一阶LPF和HPF不同设置下的频率响应。滤波器被归一化,因此最大增益为0.0dB,并且没有共振峰值。请注意,对于LPF来说,响应在奈奎斯特频率上归零,对于HPF来说,响应在直流时归零。模拟LPF的传递函数在奈奎斯特频率上没有零点,所以这代表一个错误。然而,模拟HPF的传递函数在直流处有一个零点,所以它的响应(基本上)是完全匹配的。我们将在本章末尾研究奈奎斯特频率的零点问题,并介绍一些音频滤波器的算法。
指定。
fc: 角频率
设计方程如下。
11.3.2.2 二阶 LPF 和 HPF
图11.12a和b显示了二阶LPF和HPF各种Q值设置下的频率响应。直流增益被归一化为0.0dB,导致峰值增益大于1。Q值为20时产生的峰值增益为25dB或原始增益为18;在这个增益下,信号很容易被削波。
图11.12:(a)二阶LPF和(b)HPF,fc = 500 Hz,Q = 0.707、2和20。
指定。
fc:角频率
Q:控制谐振峰值的质量系数
设计方程如下。
11.3.2.3 二阶BPF和BSF
图11.13a和b显示了二阶BPF和BSF在不同Q值设置下的频率响应。滤波器被归一化,因此最大增益为0.0dB。再次注意到,BPF的响应在Nyquist处归零--这也代表了原始模拟滤波器的一个轻微误差。
图11.13:(a)二阶BPF和(b)BSF,fc = 500 Hz,Q = 0.2、2和20。
指定。
fc:角频率
Q:控制峰值或凹槽宽度的质量系数=1/带宽
设计方程如下。
11.3.2.4 二阶巴特沃斯 LPF 和 HPF
Butterworth LPF和HPF是普通二阶HPF的专门版本。它们的Q值被固定在0.707,这是在观察到频率响应达到峰值之前它能承担的最大值。其频率响应与Q值为0.707的普通二阶滤波器基本相同。
指定。
fc: 角频率
设计方程如下