二阶Butterworth低通滤波器系数的计算

Question

带着，

采样频率:10 Freq

切断频率:1千赫

下面的差分方程的系数是如何计算的？

我知道差分方程将是这种形式，但不知道如何实际计算出系数b0，b1，b2，a1，a2的数值。

y(n)  =  b0.x(n) + b1.x(n-1) + b2.x(n-2) + a1.y(n-1) + a2.y(n-2)

我最终将在C++中实现这个LPF，但我需要知道如何首先实际计算系数，然后才能用它得到任何结果。

Answer 1

这就是你要的。ff是频率比，在您的例子中是0.1：

    const double ita =1.0/ tan(M_PI*ff);
    const double q=sqrt(2.0);
    b0 = 1.0 / (1.0 + q*ita + ita*ita);
    b1= 2*b0;
    b2= b0;
    a1 = 2.0 * (ita*ita - 1.0) * b0;
    a2 = -(1.0 - q*ita + ita*ita) * b0;

结果是：

b0=0.0674553

b1=0.134911

b2=0.0674553

a1=1.14298

a2=-0.412802

Answer 2

对于那些想知道其他答案的神奇公式从何而来的人，下面是这个例子后面的一个推导。

从Butterworth滤波器的传递函数开始

G(s) = wc^2 / (s^2 + s*sqrt(2)*wc + wc^2)

当wc是截止频率时，应用双线性z-变换，即替代s = 2/T*(1-z^-1)/(1+z^-1)。

G(z) = wc^2 / ((2/T*(1-z^-1)/(1+z^-1))^2 + (2/T*(1-z^-1)/(1+z^-1))*sqrt(2)*wc + wc^2)

T是抽样周期。

截止频率需要预先校正，以补偿z变换引入的模拟频率和数字频率之间的非线性关系：

wc = 2/T * tan(wd*T/2)

其中wd是所需的截止频率rad/s。

让C = tan(wd*T/2)，为了方便，让wc = 2/T*C。

将此替换到方程中，2/T因子就会退出：

G(z) = C^2 / ((1-z^-1)/(1+z^-1))^2 + (1-z^-1)/(1+z^-1)*sqrt(2)*C + C^2)

用(1+z^-1)^2将分子和分母相乘并展开，从而产生：

G(z) = C^2*(1 + 2*z^-1 + z^-2) / (1 + sqrt(2)*C + C^2 + 2*(C^2-1)*z^-1 + (1-sqrt(2)*C+C^2)*z^-2')

现在，将分子和分母除以常量项和分母。为了方便起见，让D = 1 + sqrt(2)*C + C^2

G(z) = C^2/D*(1 + 2*z^-1 + z^-2) / (1 + 2*(C^2-1)/D*z^-1 + (1-sqrt(2)*C+C^2)/D*z^-2')

此表格与我们正在寻找的表格相当：

G(z) = (b0 + b1*z^-1 + b2*z^-1) / (1 + a1*z^-1 +a2*z^-2)

因此，我们通过把它们等同起来，得到系数：

a0 = 1

a1 = 2*(C^2-1)/D

a2 = (1-sqrt(2)*C+C^2)/D

b0 = C^2/D

b1 = 2*b0

b2 = b0

其中，D = 1 + sqrt(2)*C + C^2，C = tan(wd*T/2)，wd是期望的截止频率rad/ s，T是采样周期s。

Answer 3

您可以使用这个链接得到具有特定采样率和频率削减的n阶Butterworth滤波器的系数。为了检验结果。您可以使用MATLAB获得系数，并与程序输出进行比较。

http://www.exstrom.com/journal/sigproc

fnorm = f_cutoff/(f_sample_rate/2); % normalized cut off freq, http://www.exstrom.com/journal/sigproc
% Low pass Butterworth filter of order N
[b1, a1] = butter(nth_order, fnorm,'low');