巴特沃斯(Butterworth)滤波器（二）_二阶巴特沃斯滤波器公式

前面讲了下巴特沃斯（Butterworth）模拟滤波器的设计，这节讲下如何将模拟滤波器离散化，变成数字滤波器，然后在单片机中实现。

我们前面讲过很多离散化方法，其中最常用的一种叫做双线性变换
s = 2 T s 1 − z − 1 1 + z − 1

$$\begin{aligned} s=\frac{2}{T_{s}}\frac{1-z^{-1}}{1+z^{-1}} \end{aligned}$$ s=Ts​2​1+z−11−z−1​​
把这个公式代入传递函数就可以得到一个Z域的差分方程。
但是如果我们直接使用双线性变换进行离散化后，得到的数字滤波器和模拟滤波器的幅频响应曲线并不一样。

可以看出数字滤波器曲线，比模拟滤波器曲线衰减的要快，如果说模拟滤波器的截止频率为10rad/s,那么数字滤波器的截止频率只有9rad/s。

这是为什么呢，因为双线性变换是近似变换，不是准确换算。那么数字滤波器的截止频率fd和模拟滤波器的截止频率fa存在着什么样的对应关系呢？

所以在设计数字滤波器的时候，先确定数字滤波器的截止频率，再通过上式计算对应模拟滤波器的截止频率，最后通过双线性变换法离散化，得到Z域的差分方程。
但值得一提的是，当fs>>fd时，即采样频率远大于截止频率时，可以得到：

模拟滤波器的截止频率等于数字滤波器的截止频率。

数字滤波器的设置步骤

1. 根据数字滤波器的截止频率计算得到对应的模拟滤波器的截止频率；
2. 根据得到模拟滤波器的截止频率，将滤波器传递函数去归一化；
3. 对模拟滤波器传函进行双线性变换；
4. 写出代码

例如：
如果我们期望的数字滤波器的截止频率为fd，对应的模拟滤波器的截止角频率为：
w a = 2 π f a = 2 f s t a n ( π f d f s )

$$\begin{aligned} w_{a}=2\pi f_{a}=2f_{s}tan(\frac{\pi f_{d}}{f_{s}}) \end{aligned}$$ wa​=2πfa​=2fs​tan(fs​πfd​​)​
去归一化只需要模拟滤波器传递函数中的s进行如下替换即可：
$$\begin{array}{r}s\to \frac{s}{w_a}\end{array}$$
最后使用双线性变换离散化，把公式代入
$$\begin{array}{r}s=2f_s\frac{1-z^{-1}}{1+z^{-1}}\end{array}$$
所以最终的变换方式就是
$$\begin{array}{r}s\to \frac{s}{w_a}=\frac{2f_s\frac{1-z^{-1}}{1+z^{-1}}}{s{f}_{s}tan(\frac{\pi f_d}{f_s})}=\frac{\frac{1-z^{-1}}{1+z^{-1}}}{tan(\frac{\pi f_d}{f_s})}\end{array}$$
举个列子：
我们以归一化的二阶巴特沃斯滤波器公式为例
$$\begin{array}{r}{H}_{(s)}=\frac{1}{s^2+1.414s+1}\end{array}$$

代入完整的变换可得：

我们令
t a n ( π f d f s ) = Ω

$$\begin{aligned} tan(\frac{\pi f_{d}}{f_{s}})=\Omega \end{aligned}$$ tan(fs​πfd​​)=Ω​

我们把上面差分方程写成代码：

const float fr = sample_freq/_cutoff_freq;
const float ohm = tanf(PI/fr);
const float c = 1.0f + 1.414f*ohm + ohm*ohm;

_b0 = ohm*ohm/c;
_b1 = 2.0f*_b0;
_b2 = _b0;

_a1 = 2.0f*(ohm*ohm - 1.0f)/c;
_a2 = (1.0f - 1.414f*ohm + ohm*ohm)/c;

for i = 1:10
	xn0 = sample(i);%采样数据
	yn0 = b0*xn0 + b1*xn1 + b2*xn2 - a1*yn1 - a2*yn2；
	yn2 = yn1；
	yn1 = yn0；
	xn2 = xn1;
    xn1 = xn0;
end
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19

这种形式有个名称叫做直接I型，我们可以将它化成框图形式

从框图中可以直观的看见使用了4个延时模块,从代码里也可以看见，我们需要保存4个过去的值；
一个二阶滤波器需要4个延时模块我们能不能化简一下呢？
我们可以把直接I型看成B(z)和A(z)两个模块串联而成的，从系统函数出发，调换其中B(z)和A(z)的顺序不影响系统输出。即：

仔细观察其实我们可以共用延时模块

这样我们只需要使用两个延时模块就行，这个叫做直接II型结构图，写成代码如下：

for i = 1:10
    d0 = sample(i) - d1*a1 -d2*a2;
    output = d0*b0 + d1*b1 + d2*b2;
    d2 = d1;
    d1 = d0;
end
1
2
3
4
5
6

这样一个巴特沃斯二阶低通数字滤波器就设置完成。

参考模型：https://download.csdn.net/download/wanrenqi/85367463
文章转载：https://zhuanlan.zhihu.com/p/357619650