跳转至

LTI 滤波器设计

可实现性、逼近准则与设计流程

滤波器对输入的不同频率分量施加不同权重，“频率”也可推广为空间频率等其他分解维度。只用过去和当前数据估计当前输出的是严格意义上的滤波器；同时使用未来数据的是平滑器，只用已有数据估计未来输出的是预测器。这里讨论的数字滤波器限定为因果 LTI 系统。

实际可实现的滤波器至少应当稳定、有限阶，并且冲激响应至少是左端有界的右边序列，才能用有限延时移成因果序列。有限阶系统函数写成

\[ H(z)=\frac{\sum_{k=0}^{M}b_kz^{-k}}{1-\sum_{k=1}^{N}a_kz^{-k}} \]

\(N=0\) 时为 FIR，否则通常为 IIR。

对有限能量因果冲激响应，课件采用的 Paley–Wiener 条件为

\[ \int_{-\pi}^{\pi}\left|\ln|H(e^{j\omega})|\right|\,\mathrm d\omega<\infty \]

因此因果稳定滤波器不能在一段非零宽度的频带上严格为零。课件还把相应结论表述为“不能在一个有限区间内恒为常数”；这里针对的是非平凡选择性响应，不包括全频带本来就恒幅的全通系统。对有限阶选择性滤波器，某段频带精确恒定、另一段又突然变为严格零的理想砖墙特性不可实现；理想低通、高通、带通和带阻都只能由因果稳定系统近似。

给定期望响应 \(H_d(e^{j\omega})\) 后，常见逼近准则有：

\[ E_2=\left[\frac{1}{2\pi}\int_{\mathcal B}|H_d(e^{j\omega})-H(e^{j\omega})|^2\,\mathrm d\omega\right]^{1/2} \]
\[ E_\infty=\max_{\omega\in\mathcal B}|H_d(e^{j\omega})-H(e^{j\omega})| \]

最大平坦准则则令 \(H\) 与 \(H_d\) 在指定频率点的若干阶导数相同。均方准则控制总体误差，最大误差准则控制最坏点，最大平坦准则把精度集中在某些关键频率附近。

设计流程包括：先给出通、阻带边缘、允许纹波和相位要求；再选择 FIR 或 IIR 形式并完成系数设计；随后作理想算术仿真；根据处理器、乘法量、存储量和有限字长效应选择实现结构；最后在实际结构和字长下重新验证。

例如模拟低通指标为 \(F_p=200\,\mathrm{kHz}\)、\(F_s=250\,\mathrm{kHz}\)，采样率 \(1\,\mathrm{MHz}\)，则数字边缘为

\[ \omega_p=0.4\pi,\qquad \omega_s=0.5\pi \]

通带峰值纹波 \(0.02\) 对应

\[ 20\log_{10}(1+0.02)\approx0.172\,\mathrm{dB},\qquad 20\log_{10}(1-0.02)\approx-0.176\,\mathrm{dB} \]

阻带峰值 \(0.01\) 对应 \(-40\,\mathrm{dB}\)。

FIR 可严格实现广义线性相位，且稳定、便于 FFT 加速，但达到同一幅度指标时阶数通常较高。IIR 常以较少参数获得更窄过渡带，但相位一般非线性，反馈也使其对系数量化和舍入更敏感。最终选择还取决于现有设计工具、成本、时延和处理平台。

广义线性相位 FIR

广义线性相位响应写成

\[ H(e^{j\omega})=A(\omega)e^{-j(\omega\alpha-\beta)} \]

其中 \(A(\omega)\) 为可正可负的实函数，群延迟恒为 \(\alpha\)。由 \(A(\omega)=A^*(\omega)\) 有

\[ H(e^{j\omega})=H^*(e^{j\omega})e^{j(2\beta-2\alpha\omega)} \]

比较 DTFT 系数得到

\[ h[n]=h[2\alpha-n]e^{j2\beta} \]

对实 \(h[n]\)，\(\beta=0\) 或 \(\pi\) 时为对称，\(\beta=\pi/2\) 或 \(3\pi/2\) 时为反对称。因而实 FIR 满足广义线性相位的常用充分条件是长度为 \(M+1\)，并满足

\[ h[n]=h[M-n] \]

或

\[ h[n]=-h[M-n] \]

此时 \(\alpha=M/2\)。注意 \(M\) 是阶数，比长度少 1。

四类实系数线性相位 FIR 为：

类型 对称性 \(M\) 奇偶 必然零点 不能实现
I 对称 偶 无固定端点零点 无额外限制
II 对称 奇 \(z=-1\) 高通
III 反对称 偶 \(z=1,-1\) 低通、高通
IV 反对称 奇 \(z=1\) 低通

I 类的频率响应可写为

\[ H(e^{j\omega})=e^{-j\omega M/2}\sum_{k=0}^{M/2}a[k]\cos(k\omega) \]

其中 \(a[0]=h[M/2]\)，\(a[k]=2h[M/2-k]\)。

II 类为

\[ H(e^{j\omega})=e^{-j\omega M/2}\sum_{k=1}^{(M+1)/2}b[k]\cos\!\left[\left(k-\frac12\right)\omega\right] \]

其中 \(b[k]=2h[(M+1)/2-k]\)。

III、IV 类把余弦换为正弦并多一个 \(j\)：

\[ H_{\mathrm{III}}(e^{j\omega})=je^{-j\omega M/2}\sum_{k=1}^{M/2}c[k]\sin(k\omega) \]
\[ H_{\mathrm{IV}}(e^{j\omega})=je^{-j\omega M/2}\sum_{k=1}^{(M+1)/2}d[k]\sin\!\left[\left(k-\frac12\right)\omega\right] \]

III 类有 \(c[k]=2h[M/2-k]\)，IV 类有 \(d[k]=2h[(M+1)/2-k]\)。

对称或反对称关系给出

\[ H(z)=\pm z^{-M}H(z^{-1}) \]

实系数又给出共轭对称。因此一般复零点 \(z_k\) 会与 \(z_k^*\)、\(z_k^{-1}\)、\((z_k^*)^{-1}\) 四点成组；位于单位圆或实轴时会退化为两点组，\(z=\pm1\) 则单独出现。这些固定零点决定了四类结构对低通、高通及 Hilbert 变换器的适用范围。

窗函数法

给定理想响应 \(H_d(e^{j\omega})\)，先取 IDTFT 得无限长 \(h_d[n]\)，再平移并加窗得到因果有限长 \(h[n]\)。在 \(h[n]\) 只允许 \(0\leq n\leq M\) 非零时，Parseval 等式表明

\[ h[n]=h_d[n]R_{M+1}[n] \]

使未加权频域均方误差最小。对截止频率 \(\omega_c\)、群延迟 \(\alpha\) 的理想低通，

\[ H_d(e^{j\omega})=\begin{cases}e^{-j\omega\alpha},&|\omega|\leq\omega_c\\0,&\omega_c<|\omega|\leq\pi\end{cases} \]
\[ h_d[n]=\frac{\sin[\omega_c(n-\alpha)]}{\pi(n-\alpha)} \]

取 \(\alpha=M/2\) 并截到 \(0\leq n\leq M\)，即可得到线性相位 FIR。

时域加窗对应频域周期卷积：

\[ H(e^{j\omega})=\frac{1}{2\pi}\int_{-\pi}^{\pi}H_d(e^{j\theta})W(e^{j(\omega-\theta)})\,\mathrm d\theta \]

窗的主瓣把理想跳变展宽成过渡带，旁瓣造成通、阻带波纹。矩形窗的 Gibbs 峰值几乎不随长度下降，\(\delta_p\approx\delta_s\approx0.0895\)，最大阻带只有约 \(-21\,\mathrm{dB}\)；增大长度主要缩窄过渡带。其经验关系为

\[ \Delta\omega\approx\frac{0.89\cdot2\pi}{M+1},\qquad \omega_c=\frac{\omega_p+\omega_s}{2} \]

长度为 \(M+1\) 的常用窗可写为：

\[ w_R[n]=1 \]
\[ w_B[n]=1-\frac{2|n-M/2|}{M} \]
\[ w_H[n]=0.5-0.5\cos\frac{2\pi n}{M} \]
\[ w_{Hm}[n]=0.54-0.46\cos\frac{2\pi n}{M} \]
\[ w_{Bl}[n]=0.42-0.5\cos\frac{2\pi n}{M}+0.08\cos\frac{4\pi n}{M} \]

它们在滤波器设计中的典型过渡宽度和最小阻带衰减为：

窗 等效 Kaiser \(\beta\) 等效过渡宽度 最小阻带衰减
矩形 0 \(1.81\pi/M\) \(21\,\mathrm{dB}\)
Bartlett 1.33 \(2.37\pi/M\) \(25\,\mathrm{dB}\)
Hann 3.86 \(5.01\pi/M\) \(44\,\mathrm{dB}\)
Hamming 4.86 \(6.27\pi/M\) \(53\,\mathrm{dB}\)
Blackman 7.04 \(9.19\pi/M\) \(74\,\mathrm{dB}\)

这些设计指标与频谱分析表中的主瓣零点宽度、最高旁瓣是不同口径，不能混用。窗频谱还常用远端旁瓣的滚降速度 \(D\,\mathrm{dB/oct}\) 描述；它与峰值旁瓣电平分别反映近端泄漏和远端泄漏，不能互相替代。

Kaiser 窗用一个额外形状参数连续调节主瓣和旁瓣：

\[ w_K[n]=\frac{I_0\!\left(\beta\sqrt{1-[(n-\alpha)/\alpha]^2}\right)}{I_0(\beta)},\qquad 0\leq n\leq M,\quad \alpha=\frac M2 \]

设目标衰减 \(A=-20\log_{10}\delta\)，经验参数为

\[ \beta=\begin{cases}0.1102(A-8.7),&A>50\\0.5842(A-21)^{0.4}+0.07886(A-21),&21\leq A\leq50\\0,&A<21\end{cases} \]
\[ M\approx\frac{A-8}{2.285\Delta\omega} \]

阶数应向上取整，并按所需 FIR 类型调整奇偶性。若 \(\delta_p=0.02\)、\(\delta_s=0.01\)、\(\omega_p=0.35\pi\)、\(\omega_s=0.45\pi\)，窗法取较严的 \(\delta=0.01\)，得到 \(A=40\,\mathrm{dB}\)、\(\beta\approx3.395\)、\(M\approx45\)、\(\alpha=22.5\)。经验式不保证一次就严格达标，仍需计算实际纹波；若高通结构类型不合适或过渡带略超限，可把长度增加 2 后重算。

频率采样与等纹波思路

频率采样法直接在 \(\omega_k=2\pi k/N\) 指定

\[ H[k]=H_d(e^{j2\pi k/N}) \]

再作 IDFT 得到 \(h[n]\)。完整频率响应由这些样本通过 Dirichlet 核插值，因此在采样点与期望值完全相等，点间误差却未直接受控。期望样本满足线性相位的相应对称条件时，得到的 FIR 也具有线性相位。理想响应突变附近的点间纹波较大，可在过渡带设置一两个经过优化的中间样本减小峰值误差。

窗法对应未加权均方逼近，频率采样法对应插值逼近；等纹波设计则求加权最大误差最小的最佳一致逼近，使各个极值点近似等高交替，通常能在同一阶数下更充分地利用允许纹波。

IIR 的模拟原型变换

模拟低通原型已有成熟的闭式或表格化设计：Butterworth 为最大平坦幅度，Chebyshev I 型在通带等纹波，Chebyshev II 型在阻带等纹波，椭圆滤波器在通、阻带都等纹波。由模拟原型设计数字 IIR 时，映射应把 \(s\) 平面虚轴映到 \(z\) 平面单位圆，并把左半平面的稳定极点映到单位圆内。

冲激响应不变法直接抽样模拟冲激响应：

\[ h[n]=T_dh_c(nT_d) \]

若

\[ H_c(s)=\sum_{k=1}^{N}\frac{A_k}{s-s_k} \]

则

\[ H(z)=\sum_{k=1}^{N}\frac{T_dA_k}{1-e^{s_kT_d}z^{-1}} \]

左半平面极点满足 \(|e^{s_kT_d}|<1\)，稳定性得以保持；频率映射 \(\omega=\Omega T_d\) 是线性的。但模拟频谱会按采样频率折叠：

\[ H(e^{j\omega})=\sum_{r=-\infty}^{\infty}H_c\!\left[j\left(\frac{\omega}{T_d}-\frac{2\pi r}{T_d}\right)\right] \]

因此它只适合高频衰减足够快的模拟原型，通常用于低通和带通，不适合高通或带阻。几何上，\(z=e^{sT_d}\) 把 \(s\) 平面每相差 \(j2\pi/T_d\) 的点映到同一个 \(z\)，这种多值映射就是混叠的复平面解释。

双线性变换为

\[ s=\frac{2}{T_d}\frac{1-z^{-1}}{1+z^{-1}},\qquad z=\frac{1+(T_d/2)s}{1-(T_d/2)s} \]

它把整条虚轴一一映到单位圆，把左半平面映到单位圆内，不产生频率混叠。

若 \(s=\sigma+j\Omega\) 且 \(\sigma<0\)，则

\[ |z|^2=\frac{(2/T_d+\sigma)^2+\Omega^2}{(2/T_d-\sigma)^2+\Omega^2}<1 \]

所以稳定模拟极点一定映到单位圆内。模拟与数字频率满足

\[ \omega=2\arctan\frac{\Omega T_d}{2},\qquad \Omega=\frac{2}{T_d}\tan\frac{\omega}{2} \]

关系非线性，设计前需对每个关键边缘预畸变：

\[ \Omega_p=\frac{2}{T_d}\tan\frac{\omega_p}{2},\qquad \Omega_s=\frac{2}{T_d}\tan\frac{\omega_s}{2} \]

双线性变换避免了混叠，但模拟域和数字域的通带、过渡带比例不再相同。模拟低通先变为模拟高通、带通或带阻再作双线性映射，与先得到数字低通原型再作相应数字频率变换的结果一致；冲激响应不变法一般没有这种可交换性。

数字频率变换与直接优化

给定数字低通原型 \(H_L(z)\)，把其中每个 \(z^{-1}\) 替换为稳定全通映射 \(G(Z^{-1})\)，可得

\[ H_d(Z)=H_L(z)\big|_{z^{-1}=G(Z^{-1})} \]

映射需把单位圆映到单位圆、把单位圆内外关系保持不变，并且是便于实现的有理函数。一般形式可写为

\[ G(Z^{-1})=\pm\prod_{k=1}^{K}\frac{Z^{-1}-\alpha_k}{1-\alpha_kZ^{-1}},\qquad |\alpha_k|<1 \]

若改用课件的正幂变量记号，并令 \(D(Z)=\prod_k(Z-\alpha_k)\)，同一关系可写为

\[ G(Z)=\pm\frac{D(Z)}{Z^KD(Z^{-1})} \]

分子、分母系数互为反序，显式呈现全通结构。

一阶低通到低通变换为

\[ G(Z^{-1})=\frac{Z^{-1}-\alpha_1}{1-\alpha_1Z^{-1}},\qquad \alpha_1=\frac{\sin[(\theta_p-\omega_p)/2]}{\sin[(\theta_p+\omega_p)/2]} \]

\(\theta_p\) 是原型边缘，\(\omega_p\) 是目标边缘。低通到高通变换为

\[ G(Z^{-1})=-\frac{Z^{-1}+\alpha_2}{1+\alpha_2Z^{-1}},\qquad \alpha_2=-\frac{\cos[(\theta_p+\omega_p)/2]}{\cos[(\theta_p-\omega_p)/2]} \]

更高阶全通映射可生成带通、带阻、多通带和多阻带响应。

对不规则幅频指标，也可直接把系统写成二阶节级联

\[ H(z)=A\prod_{k=1}^{K}\frac{1+a_kz^{-1}+b_kz^{-2}}{1+c_kz^{-1}+d_kz^{-2}}=A\,G(z) \]

在一组频率 \(\omega_i\) 上最小化幅度误差平方和：

\[ E=\sum_{i=1}^{P}\left[A|G(e^{j\omega_i})|-|H_d(e^{j\omega_i})|\right]^2 \]

固定各二阶节参数后，最优总增益为

\[ A_{\mathrm{opt}}=\frac{\sum_i|G(e^{j\omega_i})|\,|H_d(e^{j\omega_i})|}{\sum_i|G(e^{j\omega_i})|^2} \]

其余 \(4K\) 个参数满足一组非线性方程，需要迭代求解。优化结果不自动保证极点在单位圆内；可把单位圆外极点镜像入内，并乘相应全通因子保持幅度。这类方法适合没有标准原型的不规则响应。

实现结构

信号流图由加法、常数乘法和单位延迟组成。对单输入单输出流图，反转所有支路方向、保持支路系数不变并交换输入与输出，所得转置流图具有相同系统函数。转置形式常把一个多输入加法器改成多个两输入加法器，也会改变内部节点动态范围和有限字长误差。

IIR 直接 I 型把分子 FIR 和分母递归部分分别实现，需要约 \(M+N\) 个延迟；交换两部分次序不影响总响应。直接 II 型合并两组延迟，内部变量满足

\[ w[n]=x[n]+\sum_{k=1}^{N}a_kw[n-k] \]
\[ y[n]=\sum_{k=0}^{M}b_kw[n-k] \]

只需 \(\max(M,N)\) 个延迟，是规范型实现，但内部状态的动态范围与噪声增益可能更大。

把实系数有理函数分成一、二阶节可得级联形式：

\[ H(z)=\prod_{k=1}^{N_s}\frac{b_{0k}+b_{1k}z^{-1}+b_{2k}z^{-2}}{1-a_{1k}z^{-1}-a_{2k}z^{-2}} \]

它便于独立调整零极点、模块化实现，并可通过零极点配对和节次排序控制溢出与舍入噪声。部分分式展开得到并联形式，各支路是一、二阶节；支路误差不再沿后续节级联放大，也便于并行和模块化实现。

FIR 直接型就是抽头延迟线：

\[ y[n]=\sum_{k=0}^{M}h[k]x[n-k] \]

转置型由流图转置得到。若 \(h[M-n]=\pm h[n]\)，先把对称位置的输入相加或相减，再乘一次共同系数，可把乘法器数约减半。FIR 也可按零点分为标准线性相位节：

\[ H_1(z)=1\pm z^{-1} \]
\[ H_2(z)=1-2\cos\theta\,z^{-1}+z^{-2} \]
\[ H_3(z)=1-\left(r+\frac1r\right)z^{-1}+z^{-2} \]

以及由共轭镜像四点组形成的四阶回文节

\[ H_4(z)=1+bz^{-1}+cz^{-2}+bz^{-3}+z^{-4} \]

其中

\[ b=-2\left(r+\frac1r\right)\cos\theta,\qquad c=r^2+r^{-2}+4\cos^2\theta \]

频率采样结构从 \(H[k]=\operatorname{DFT}\{h[n]\}\) 出发：

\[ H(z)=\frac{1-z^{-N}}{N}\sum_{k=0}^{N-1}\frac{H[k]}{1-W_N^{-k}z^{-1}} \]

窄带滤波器只有少量 \(H[k]\) 非零时，所需支路很少，且各系数直接对应频率样本，适合模块化和时分复用。但该表达式含单位圆上的递归极点，理论上由前置零点精确相消，有限字长下相消可能失配；可把半径改为略小于 1 的 \(r\)：

\[ H(z)=\frac{1-r^Nz^{-N}}{N}\sum_{k=0}^{N-1}\frac{H[k]}{1-W_N^{-k}rz^{-1}} \]

以保证递归支路稳定。结构选择不能只看输入输出传函；延迟数、并行能力、系数量化灵敏度、节点动态范围、舍入噪声和极限环都可能不同。

评论