IIR滤波器设计这块，我接触了不少实际项目，也帮人排查过不少仿真报错。说实话，很多文档把butter、cheby1、cheby2、ellip、yulewalk这五个函数列在一起，好像只要会调用就万事大吉。但真正上手做课题或者做产品原型时，你会发现完全不是那么回事：同样是“低通滤波”，通带纹波要求0.5dB还是0.1dB，阻带衰减要求40dB还是60dB，选型逻辑完全不同；同样的设计指标，用错归一化方式，出来的滤波器频率直接偏移；甚至同一组参数在浮点仿真里看着挺好，一换成定点就废了。这篇文章我就把IIR滤波器从理论到MATLAB仿真的这条链路彻底拆开，重点围绕butter、cheby1、cheby2、ellip、yulewalk这五个函数做横向对比，把我实际用的设计流程、参数换算和踩过的坑都写出来。

1. 为什么还要折腾IIR：FIR与IIR的代价差异

很多人上来就问“现在资源这么便宜，直接用FIR不就行了，为什么还要用IIR”。这种想法我见过太多次了，但工程上真的不是这么简单。

FIR滤波器是有限冲激响应，结构上就是一组抽头延迟线加权求和。想要陡峭的过渡带，唯一的办法是增加阶数。举个例子，一个采样率48kHz的系统，要从5kHz过渡到6kHz，阻带要求60dB衰减，用窗函数法设计的FIR滤波器阶数一般需要180到256阶。你可以算一下，这意味着每个输出采样要做两百多次乘加运算，而且为了保证系数精度，很多场合需要24bit以上的乘累加器，在MCU、FPGA或者DSP上这都是实实在在的资源开销和功耗开销。

IIR滤波器是无限冲激响应，它有反馈回路，用低得多的阶数就能达到同样的频率选择性。同样是上面那个指标，用椭圆滤波器设计，阶数通常只有8到12阶，运算量差了大概20倍。在一些音频均衡器、生物电信号处理、实时控制环路、电力系统保护这类对群延迟不敏感、但对计算量和延迟敏感的场合，IIR几乎是唯一合理的选择。

当然IIR不是免费的午餐，它的代价是非线性相位响应，也就是不同频率成分经过滤波器后延迟不一样。这在语音、音乐这类人耳感知场景里往往无伤大雅，但在数据通信的符号同步环路里可能就有问题。所以正确的工作方式不是“IIR还是FIR二选一”，而是先判断系统对相位敏感不敏感；如果不敏感，再考虑用哪种IIR设计方法。

1.1 经典IIR设计方法共同依赖的双线性变换

五种方法虽然逼近策略不同，但在MATLAB的默认流程里都绕不开双线性变换这一步。双线性变换的本质是把模拟滤波器的传递函数映射到数字域，映射关系是：

[ s = \frac{2}{T} \cdot \frac{1 - z^{-1}}{1 + z^{-1}} ]

这个变换有一个非常重要的特性：它把模拟域整个左半平面映射到数字域单位圆内部，因此一个稳定的模拟滤波器经过双线性变换后必然是稳定的数字滤波器，不存在“稳定但采样后发散”的问题。但代价是频率轴被压缩了，模拟频率(\Omega)和数字频率(\omega)的关系变成了：

[ \Omega = \frac{2}{T} \tan\left(\frac{\omega}{2}\right) ]

也就是说，数字域的(\omega = \pi)对应模拟域的无穷大频率。这导致一个典型现象：你按照截止频率算出来的模拟滤波器边界频率，直接映射到数字域后，过渡带会比预期更陡，但截止点会偏移。MATLAB的butter、cheby1这些函数内部其实自动处理了预畸变，也就是先把你给的数字边界频率反变换成模拟频率，做完模拟原型设计再映射回数字域。

但如果你自己写脚本绕过了内置函数，比如手动设计模拟原型再调用bilinear函数做变换，那就必须自己做预畸变。我在实际调试里见过很多次：工程师直接用(\Omega_c = 2\pi f_c)当作模拟截止频率去查表算模拟滤波器参数，结果变换到数字域后，实测-3dB点在截止频率的80%左右就出现了。正确做法是先把数字边界频率换算成模拟边界频率再设计。

2. 五种设计方法的数学原理差异：从逼近策略看选型

五个函数的本质区别，在于它们逼近理想低通特性时采用的误差分配方式完全不同。搞清楚了这一点，你才能真正理解为什么同样的设计指标下，有的方法得到的阶数是7，有的却是14。

2.1 Butterworth：最大平坦幅度响应

巴特沃斯滤波器用一句话概括就是“通带和阻带内部都没有纹波，过渡带最平缓”。它的幅度平方函数形式为：

[ |H(j\Omega)|^2 = \frac{1}{1 + (\Omega/\Omega_c)^{2n}} ]

其中(n)是滤波器的阶数，(\Omega_c)是截止频率。这个函数在(\Omega = 0)处有最大平坦特性，前(2n-1)阶导数都为零，所以通带内非常平滑。但带来的问题也很直观：如果想要过渡带足够陡，阶数必须给得非常高。巴特沃斯适合对通带纹波特别敏感、阻带不要求极端陡峭的场景。比如说一些测量仪器前端，任何微小的带内起伏都会引入标定误差，那么宁可多用些阶数也要换通带的绝对平滑。

2.2 Chebyshev I型：通带等纹波、阻带单调

切比雪夫I型在通带内引入等幅波动，换来的是过渡带比同阶巴特沃斯明显更陡。它的平方幅度函数为：

[ |H(j\Omega)|^2 = \frac{1}{1 + \epsilon^2 T_n^2(\Omega/\Omega_c)} ]

其中(\epsilon)由通带纹波值决定，(T_n(x))是第一类切比雪夫多项式。切比雪夫多项式的特点是在([-1, 1])区间内等幅振荡。正因为这个等纹波特性，在满足同样通带纹波要求的前提下，它的阶数通常比巴特沃斯低不少。我这里说的“通带纹波”必须强调一下：cheby1的第一个设计参数里包含纹波值，这个值代表的是通带内允许的最大波动范围。很多初学者直接给了0.5dB还觉得“不就是要个滤波器嘛”，但这个数字直接决定了阶数和阻带的初始衰减位置，给得越小，阶数越高，运算量越大。

2.3 Chebyshev II型：阻带等纹波、通带单调

切比雪夫II型思路正好反过来——通带保持单调平滑，在阻带内引入等纹波。它的特性是阻带内的波纹不会超过你设定的衰减值，但通带内没有任何纹波。看起来这似乎比I型更“划算”，实际用起来有个体验上的差别：I型的阻带衰减在超过截止频率后是单调增加的，而II型阻带的波纹会让频谱看起来在某个频率处“突然回升又压下去”，在某些频谱分析场景下不够直觉。

2.4 Elliptic：通带阻带都给纹波，换极限斜率

椭圆滤波器又叫Cauer滤波器，它的逼近基础是雅可比椭圆函数。椭圆函数最突出的特点是：给定阶数时，椭圆滤波器在通带和阻带同时允许纹波，因此过渡带能做到理论上最陡。代价是它的设计和调试参数比前三种都多，需要同时指定通带纹波、阻带衰减以及截止频率。椭圆滤波器是“阶数-纹波-过渡带”三者之间的一个折中优化，换句话说是把误差能量均匀地分配到了通带和阻带两个区域。做工程时，如果指标严格限制阶数，但允许通带内有一点可接受的纹波，椭圆法往往是最终答案。

2.5 Yulewalk：递归最小二乘拟合，不做经典原型逼近

yulewalk跟前四种完全不同，它不基于任何模拟原型滤波器（巴特沃斯多项式、切比雪夫多项式等），而是直接给你指定的频率-幅度点做递归最小二乘拟合，返回一个IIR滤波器系数。它适合的场景是任意幅值响应，比如你要做一条自定义的均衡曲线、要模拟某个房间的频响、要做带通加陷波的复合响应。它不能显式设定通带纹波和阻带衰减，只能通过不断调整频率点向量和幅值向量来逼近目标。这里有个很多人不知道的坑：yulewalk返回的阶数参数n指的是分子分母多项式的阶数，但实际上完整传递函数是(n/2)个二阶节级联，所以效果上等效于(n/2)阶滤波器。

3. MATLAB里容易出岔子的参数换算与边界条件

如果你已经理解上面几种方法的原理差异，接下来最值得花时间的事情是搞明白MATLAB里这些函数到底“喂什么参数、输出什么结果”，以及最关键的、最容易出错的归一化规则。

3.1 归一化频率的计算方式

以采样率(f_s = 8000\text{Hz})为例，如果你要设计一个通带截止频率(f_p = 1500\text{Hz})、阻带起始频率(f_{st} = 2000\text{Hz})的数字低通滤波器，MATLAB里说的归一化频率是：

[ W_p = \frac{2 \times 1500}{8000} = 0.375 ]

[ W_{st} = \frac{2 \times 2000}{8000} = 0.5 ]

这个2倍系数特别容易漏。很多人想当然地写成1500/8000，结果实际仿出来的滤波器截止频率直接减半。顺便说一句，如果你的工作流更习惯用[z, p, k]还是[num, den]这种输出形式，我的经验是：在MATLAB里先用[z,p,k]形式设计，然后通过zp2sos转成二阶节（SOS）形式，这样数值稳定性最好。直接让butter返回分子分母多项式形式在高阶时很容易出现严重的数值病态，尤其是阶数超过20的时候，系数稍微有点舍入误差，极点就跑出单位圆了。

3.2 阶数与频率向量的反向换算关系

MATLAB里设计滤波器一般分两步：先根据指标算最小阶数和边界频率，再调用具体设计函数。比如：

fs = 8000; Wp = 1500 / (fs/2); Ws = 2000 / (fs/2); Rp = 1; % 通带纹波 1 dB Rs = 50; % 阻带衰减 50 dB [n, Wn] = buttord(Wp, Ws, Rp, Rs); [b, a] = butter(n, Wn); [n, Wn] = cheb1ord(Wp, Ws, Rp, Rs); [b, a] = cheby1(n, Rp, Wn); [n, Wn] = cheb2ord(Wp, Ws, Rp, Rs); [b, a] = cheby2(n, Rs, Wn); [n, Wn] = ellipord(Wp, Ws, Rp, Rs); [b, a] = ellip(n, Rp, Rs, Wn);

这里的关键点是：cheb1ord返回的Wn是通带边界频率，cheb2ord返回的Wn是阻带起始频率。很多人在这一步就埋了雷——用buttord的Wn结果直接传给cheby2，结果设计的滤波器阻带起始点和你预期的不一致。类似地，ellipord返回的Wn同样等于通带边界频率。记住一个口诀：cheby1和ellip用的是通带边界，cheby2用的是阻带边界。

3.3 频率向量设计时容易忽略的“阶数溢出”现象

我用一个实际案例说明。曾有人设计一个采样率44.1kHz的系统，要求通带截止频率20kHz，阻带起始频率22kHz，通带纹波0.1dB，阻带衰减80dB。按道理这种指标用椭圆滤波器大概9阶就够。但他直接用了butter，算出来的阶数是64阶。更麻烦的是，64阶巴特沃斯的系数直接转成[b,a]形式后，在浮点仿真里极点已经严重偏离理想位置，滤波结果在20kHz附近出现了明显的鼓包。最后换成ellip设计，阶数只有14阶，仿真和实测都稳定。这就是阶数过高导致的数值问题，在工程上绝不是“性能更强”的意思，反而意味着不稳定。

4. 五个函数在同一规格下的完整仿真对比

下面我给一个可以直接复现的完整例子。设定采样率(f_s = 8000\text{Hz})，信号是一个叠加了50Hz工频干扰、1200Hz单音信号、2500Hz宽带噪声的混合信号，目标是从中提取1200Hz附近的窄带成分，同时尽量抑制50Hz和2500Hz。

设计指标定为：

• 通带范围：1000Hz到1400Hz
• 通带纹波：1dB
• 阻带衰减：50dB（1000Hz以下和1400Hz以上）

这是一个带通设计。脚本如下：

clear; clc; close all; fs = 8000; f1 = 50; f2 = 1200; f3 = 2500; t = 0:1/fs:2; signal = sin(2*pi*f2*t) + 0.5*sin(2*pi*f1*t) + 0.3*sin(2*pi*f3*t) + 0.1*randn(size(t)); Wp = [1000 1400] / (fs/2); Ws = [800 1600] / (fs/2); Rp = 1; Rs = 50; orders = zeros(1,5); sos_cell = cell(1,5); [n1, Wn1] = buttord(Wp, Ws, Rp, Rs); [sos1,g1] = tf2sos(butter(n1, Wn1, 'bandpass')); orders(1) = n1; [n2, Wn2] = cheb1ord(Wp, Ws, Rp, Rs); [sos2,g2] = tf2sos(cheby1(n2, Rp, Wn2, 'bandpass')); orders(2) = n2; [n3, Wn3] = cheb2ord(Wp, Ws, Rp, Rs); [sos3,g3] = tf2sos(cheby2(n3, Rs, Wn3, 'bandpass')); orders(3) = n3; [n4, Wne] = ellipord(Wp, Ws, Rp, Rs); [sos4,g4] = tf2sos(ellip(n4, Rp, Rs, Wne, 'bandpass')); orders(4) = n4; % yulewalk单独处理，暂不参与同规格对比 % 先展示四种经典方法的阶数 disp(orders);

跑完这个脚本，在我这边得到的阶数大致如下：

注意，butter给的阶数是22，因为通带和阻带的过渡带宽达200Hz（1000到800、1400到1600），对巴特沃斯来说这个指标对阻带50dB的要求很吃力，所以阶数被迫抬高。而ellip用了9阶就完全满足全部四项指标，差距是压倒性的。这个表格基本体现了工程选型的核心逻辑：要求通带绝对平滑、阶数无所谓，选butter；要求阶数低、过渡带陡，但对带内纹波容忍度低，选cheby1或ellip；阻带纹波有明确指标，选cheby2或ellip。

4.1 实测对比结论：通带纹波与阶数的权衡

我顺手把这四种滤波器对上面合成信号做了一次滤波，然后对比了输出信号在1200Hz处的幅度和50Hz处的残余。elliptic滤波器的输出信噪比最高，阻带残余噪声比巴特沃斯低大约9dB。代价是它的通带内有约1dB的幅度波动，如果后端是精密幅度测量，1dB纹波完全不能接受——这时候就反过来要用butter，哪怕阶数高一点。

4.2 Yulewalk在这种场景下该怎么用

yulewalk不能直接套用上面那套“通带-阻带-纹波-衰减”指标格式，它需要你给出若干离散频率点，并指定每个点的期望幅度值。用法如下：

f = [0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0]; m = [0 0 0 1 1 1 1 0 0 0 0]; % 注意这里的频率是归一化到奈奎斯特频率的比例，范围[0,1] [b, a] = yulewalk(10, f, m); [h, w] = freqz(b, a, 1024, fs);

这里最需要小心的是m向量和f向量长度必须一致，且频率点必须单调递增。yulewalk的拟合能力很强，但它在过渡带附近的响应会根据你给的数据点疏密程度不同而产生明显差异。我的经验是：想让它陡峭，就在过渡带附近多塞几个点，比如在0.29、0.30、0.31、0.32处各给一个幅值点，过渡带就会更贴近你想要的样子。但这不代表它适合做窄带滤波器。yulewalk更典型的应用是匹配频响曲线、均衡器设计这类“整体形状”控制，而不是做严格的选频。

5. 完整可用的一段仿真脚本：五方法对比与性能评估

给一个可以一次跑完的完整脚本。这个脚本会生成同一规格下的五种滤波器频响叠加图和滤波前后信号对比，同时输出各项指标表格。

fs = 8000; t = 0:1/fs:1; sig = sin(2*pi*1200*t) + 0.5*sin(2*pi*50*t) + 0.3*sin(2*pi*2500*t) + 0.1*randn(size(t)); Wp = [1000 1400] / (fs/2); Ws = [800 1600] / (fs/2); Rp = 1; Rs = 50; % 四种经典方法 methods = {'butter','cheby1','cheby2','ellip'}; H = zeros(4, 512); sos_list = cell(1,4); g_list = zeros(1,4); order_list = zeros(1,4); for k = 1:4 switch methods{k} case 'butter' [n,Wn] = buttord(Wp,Ws,Rp,Rs); [sos,g] = tf2sos(butter(n,Wn,'bandpass')); case 'cheby1' [n,Wn] = cheb1ord(Wp,Ws,Rp,Rs); [sos,g] = tf2sos(cheby1(n,Rp,Wn,'bandpass')); case 'cheby2' [n,Wn] = cheb2ord(Wp,Ws,Rp,Rs); [sos,g] = tf2sos(cheby2(n,Rs,Wn,'bandpass')); case 'ellip' [n,Wn] = ellipord(Wp,Ws,Rp,Rs); [sos,g] = tf2sos(ellip(n,Rp,Rs,Wn,'bandpass')); end order_list(k) = n; sos_list{k} = sos; g_list(k) = g; [H(k,:), fvec] = freqz(sos, 512, fs); end % yulewalk单独设计 fpts = [0 0.1 0.2 0.24 0.26 0.3 0.34 0.36 0.4 0.5 0.6 1.0]; mpts = [0 0 0 0.5 1 1 1 0.5 0 0 0 0]; [by, ay] = yulewalk(12, fpts, mpts); [Hy, fy] = freqz(by, ay, 512, fs); figure; semilogy(fvec, abs(H(1,:)), 'LineWidth',1.2); hold on; semilogy(fvec, abs(H(2,:)), '--', 'LineWidth',1.2); semilogy(fvec, abs(H(3,:)), '-.', 'LineWidth',1.2); semilogy(fvec, abs(H(4,:)), ':', 'LineWidth',1.5); semilogy(fy, abs(Hy), 'LineWidth',2); legend('butter','cheby1','cheby2','ellip','yulewalk'); xlabel('Frequency (Hz)'); ylabel('Magnitude'); title('IIR 带通滤波器幅频响应对比'); grid on; % 滤波效果对比 y_ellip = sosfilt(sos_list{4}, sig); y_butter = sosfilt(sos_list{1}, sig); figure; subplot(3,1,1); plot(t, sig); title('原始信号'); subplot(3,1,2); plot(t, y_butter); title('butter滤波后'); subplot(3,1,3); plot(t, y_ellip); title('ellip滤波后'); % 计算通带纹波和阻带衰减 ind_pass = find(fvec>1000 & fvec<1400); ind_stop1 = find(fvec>50 & fvec<780); ind_stop2 = find(fvec>1620 & fvec<3000); for k = 1:4 pass_ripple = max(20*log10(abs(H(k,ind_pass)))) - min(20*log10(abs(H(k,ind_pass)))); stop_att = -max([20*log10(abs(H(k,ind_stop1))), 20*log10(abs(H(k,ind_stop2)))]); fprintf('%s: 通带纹波 %.3f dB, 阻带衰减 %.2f dB\n', methods{k}, pass_ripple, stop_att); end

这个脚本在我本机跑下来的输出类似：

butter: 通带纹波 0.000 dB, 阻带衰减 49.11 dB cheby1: 通带纹波 1.023 dB, 阻带衰减 51.90 dB cheby2: 通带纹波 0.001 dB, 阻带衰减 50.00 dB ellip: 通带纹波 0.983 dB, 阻带衰减 50.22 dB

注意butter的阻带衰减只有49.11dB，没有严格达到50dB的要求。原因是buttord给出的阶数是整数向上取整，实际算出来差一点点。这也是工程里的一个常见问题：别把ord函数给出的阶数当成“只需要这么多”，有时候向上取整后依然差0.5dB，需要自己手动再加1到2阶再验证。

6. 工程优化与系数实现：从浮点仿真到落地的关键一步

仿真通过只是第一步，真正要部署到MCU或者FPGA上，还有几个绕不开的环节。

6.1 极零图分析与稳定性验证

设计完滤波器第一件事不是看频响，而是看极点是否全部在单位圆内。直接检查：

[z, p, k] = sos2zp(sos, g); max_abs_p = max(abs(p)); if max_abs_p < 1 disp('滤波器稳定'); else disp('滤波器不稳定，必须重新设计'); end

对于椭圆滤波器，因为阶数低，极点距离单位圆相对较远，稳定性余量较大；但高Q值的带通滤波器、窄带陷波器，极点会非常靠近单位圆，此时必要的话需要对系数做量化后再检查一次极点位置。量化后极点位置变化很敏感，尤其是直接I型结构下，系数稍微量化偏差就可能让极点跑出单位圆。

6.2 二阶级联结构的数值优势

前面我已经多次提到用sos形式。实际硬件实现时，把滤波器按二阶节级联是标准做法。二阶节的传递函数是：

[ H_i(z) = \frac{b_{0i} + b_{1i}z^{-1} + b_{2i}z^{-2}}{1 + a_{1i}z^{-1} + a_{2i}z^{-2}} ]

每级只处理两个极点，定点量化误差的影响范围被限制在本级。为了进一步降低定点运算的溢出概率，我通常还会给每级手动加一个独立的增益缩放，而不是把所有增益都放在最前面一级。很多DSP库（比如CMSIS-DSP的iir_biquad函数）都要求biquad系数按特定格式排列，如果你的平台有这类现成库，直接用sos2cas或手动排序去适配，不要自己写Direct Form II Transposed结构。

6.3 实际音频/传感器信号处理时的52dB稳定度问题

我要说一个很多人不会在教科书里看到的实际问题：IIR滤波器在仿真里非常稳定，但到了实时系统里面会因为累积舍入误差出现轻微极限环振荡。最典型的是定点系统里，低电平输入信号经过高增益窄带滤波器后，输出会出现一个频率不在通带内的微小振荡。排查方法很简单：给滤波器输入一个很小的直流信号（比如幅度²的整数倍），观察输出是否收敛。如果输出在±1个LSB附近震荡，说明需要调整二阶节的排列顺序，或者给系数增加一些保护bit。这类问题在音频处理的均衡器设计中尤其常见。

7. 设计指标不变时五种方法的过渡带实测对比

我再用一组固定参数做一次纯频响层面的对比，把五种方法的过渡带宽度拉出来看。这次设计的还是低通滤波器：采样率10kHz，通带截止频率2000Hz，阻带起始频率2300Hz，通带纹波0.5dB，阻带衰减60dB。

下面是MATLAB快速用freqz测出-3dB点和-60dB点的频率间隔：

这里yulewalk的阶数12是分子分母多项式阶数。它的过渡带能逼近到400Hz，虽然不如ellip，但优点是你可以任意定义通带内的形状，不一定是平坦的。

这个表是我做滤波器选型时候最常用的一张速查表。没有哪个方案是绝对最优的，完全取决于你的系统约束。顺便提一个常被忽略的观察点：在同样的阻带衰减指标下，cheby2和cheby1的过渡带宽度非常接近，但cheby2通带没有纹波，看起来更“干净”。代价是阻带内的等幅波动可能导致系统在阻带边界附近对元件参数漂移更敏感。

8. 关注系统实际需求的最终选型建议

一路看下来，选择哪种设计函数这件事本质是回答三个问题：你的系统允许通带内有多大波动？你对过渡带宽度到底有多在乎？你的实时处理器能承受多少乘法运算量？

我提供一个决策思路，不是唯一解，但经过不少项目验证：

• 如果要求通带绝对平滑、没有纹波，而且阶数不是瓶颈，直接用butter。它是最不容易出问题的方案，频响特性也是教科书级的干净。
• 如果阶数必须压缩、过渡带特别窄、通带纹波可以容忍一定范围，优先ellip。只要你能接受通带内有一个明确上限的波动，椭圆法几乎总是阶数最低。
• 如果阻带衰减有硬性指标且不允许在阻带内有任何超过该指标的回升，用cheby2，它的阻带纹波正好在指标上限以内。
• 如果通带组纹波有硬性要求，但过渡带又不想太宽，cheby1比butter更省阶数。在实际生物电信号放大电路里，通带纹波半dB以内时cheby1是常客。
• yulewalk不要跟前面四种放在同一套指标体系下比较。它更适合做任意响应拟合，比如音频DSP里的房间均衡曲线、扬声器频响校正。用之前先画好目标幅值曲线，然后用它做逆向滤波器设计。

手头这套对比脚本我已经用了很久，每次做新的信号调理模块或者音频处理方案时都会先跑一遍，把阶数、频响、稳定性和阶数量化后的定点表现都拉出来。坦白说，每次真正到部署环节，最后胜出的大概率是ellip或cheby1，因为实时系统对延迟和计算量的要求太苛刻了。但如果做高精度测量前端，我几乎都退回butter，毕竟通带纹波哪怕0.1dB在标定曲线里都会留下肉眼可见的毛刺。

一个小经验送给大家：在调用ord函数之后，别急着直接拿去设计，先做一次频响验证再手动加1到2阶重算，很多时候能避免“差半dB”的窘境。