matlab 滤波器设计 coe_电机控制：滤波篇之巴特沃兹滤波器(BF) (2)

目录

一，巴特沃兹滤波器介绍
二，低通滤波中细节问题详解
三，巴特沃兹滤波器的s域求解
四，巴特沃兹滤波器离散化
五，Matlab仿真验证
六，结论

一、巴特沃兹滤波器(Butterworth Filter, BF)介绍

在电机控制中，为了进一步抑制高频噪声，科研人员也会采用一些高阶低通滤波器来对电流采样信号的高频噪声进行抑制，常用的一种滤波器为：巴特沃兹滤波器。除了这种滤波器，也存在如贝塞尔、切比雪夫滤波器等。

这里我先引用下文中TI的参考指南：

https://www.ti.com/cn/lit/an/zhct242/zhct242.pdf?ts=1605165773236&ref_url=https%253A%252F%252Fwww.google.com%252F​www.ti.com

为了使高频信号衰减，所有数据采集系统都在模数转换器（ADC）前面有一个抗锯齿（低通）滤波器，或者在数模转换器（DAC）后面有一个抗镜像 （低通）滤波器。那么如何选取合适的滤波器呢？这就是细节化的功夫活了。

我们需要了解这三种低通滤波器的幅频和相频特性，虽然这三种滤波器都是低通滤波，但是在幅频和相频特性上还是显著不同：

来源http://mcu.eetrend.com/content/2018/100012078.html

来源http://mcu.eetrend.com/content/2018/100012078.html

从上述图中，理想地，我们希望黑色虚线可以作为低通滤波器，但实际中，是不存在这种ideal的。因此退而求其次，可以发现贝塞尔滤波器具有线性相位特性，便于补偿，然而在幅频特性中阻带下降最慢。切比雪夫滤波器在幅频特性中阻带下降最快，然而在通带含有纹波，巴特沃兹滤波器在通带最平坦，阻带速度相对适中。根据上述幅频相频特性，进一步总结如下特性：

TI 讲义

为了了解滤波器具有特性的本质，这里面，仅对巴特沃兹滤波器进行介绍。

巴特沃兹滤波器的增益

:

(1)

其中，

为传递函数增益；

为输入信号的频率；

为截至频率即-3db的频率；

为滤波器阶数。当

，一阶的巴特沃兹滤波器也就是我们常用的类似于RC的低通滤波器。

二、低通滤波中细节问题详解

• 问题一：当
  ,滤波器为一阶低通滤波器，-3db为什么对应系统的截止频率？

在一阶巴特沃兹滤波器中，其传递函数可以表示如下：

(2)

当

, 以下成立：

（3）

根据上述公式，当信号频率到达截止频率的时候，幅频特性被衰减到0.707，也就是-3db。

也可以理解为系统闭环带宽

：结论参考我下述文献中的

理论分析法2测量带宽： K.0k：理解电机控制系统中的带宽问题​zhuanlan.zhihu.com

• 问题二：系统上升时间
  与系统带宽
  的关系式的由来？

我们都知道下述公式成立：

（4）

其中，上升时间为

，

为系统的带宽。

为了对其进行分析，首先，写一个含有带宽的低通滤波器通式：

(5)

我们对Eq.（5）进行s域到时域转换（这里也可以直接进行s域反变换，即另一种方法）：

(6)

其中，

分别表示系统是输出和输入信号。

我们对上述线性非齐次常系数微分方程Eq. (6)进行求解：

(7)

其中，

为系统初值，这里默认为零,

为阶跃输入信号；系统上升时间

定义为系统输出信号到达输入信号的10%-90%：

假设

,幅值在输入信号的10%-90%的输出信号可以表示为：

(8)

其中，幅值在输入信号的10%-90%的输出信号所对应的时间为

,即：

(9）

同时上升时间

可以表示为：

(10)

对上述Eq. (8)两个公式做除法：

(11)

做ln (Matlab 中的log变换):

(12)

根据上述Eq. (12)，利用matlab求解：

>> k=log(9)/(2*pi)

k = 0.3497

我们获得，

,即

（13）

值得注意的是，这个关系式仅在一阶微分方程适用，高阶微分方程在0.4附近。

三、巴特沃兹滤波器的s域求解

巴特沃兹滤波器中拥有在通带处具有最大平坦幅度，且截至频率后陡峭，这一特性来自于增益

的合理选取:

(1)

我对滤波器设计的理解是由于幅频特性是约束通带passband，阻带stopband最好的描述，因此不需要从s、z域进行入手设计，而是通过在幅频特性上对其进行设计，并转化成s域或者z域。

任意一个无零点型的传递函数都可以写成上述Eq. (1)形式，我们可以取一部分式子如：

（14）

当阶数

越大时，

变化越大，即下图所示：

因此当我们确定了阶数和截至频率，即系统增益

后，我们可以确定滤波器的增益

。然后可以进行转化成s域，假设我们需要构造一个二阶通式：

(15)

求上述式子幅值：

(16)

同时根据Eq. (1)，改写为：

(17)

根据Eqs. (16-17)，可以计算出参数

:

(18)

我们需要创造一个稳定系统，需要满足系统的所有极点位于左半轴，即:

(19)

代入Eq. (19)于Eq. (15)，我们构造的二阶s域巴特沃兹滤波器传递函数如下所示：

(20)

例如：我们需要搭建一个二阶巴特沃兹滤波器，且截至频率

为100Hz，增益

，根据Eq. (20), 我们可以获得二阶巴特沃兹滤波器的传递函数为：

(21)

同理，在求解3-4阶或者更高阶，也可以采用这种方法进行构造。

四、巴特沃兹滤波器离散化

由于实际系统为数字化系统，因此在构建了连续域后，我们需要进行离散化，这里我采用Tustin变换（也称为双线性变换）先进行到z域，然后在z域在转变成离散域，即：

(22)

可参考这篇文章：

K.0k：PI电流环在数字化实际系统中的完整解析​zhuanlan.zhihu.com

例如，我们对Eq. (21)进行离散化：

首先，代入Eq. (2)于Eq. (21)，计算出z域的传递函数：

（23）

然后对Eq. (23)进行离散域变换（Matlab代码如下）：

function filter_out     = filter(filter_in, Tk)
persistent Y_K1 Y_K2 Y_K3 X_K1 X_K2 X_K3;
if isempty(X_K3), X_K3 = 0; end
if isempty(X_K2), X_K2 = 0; end
if isempty(X_K1), X_K1 = 0; end
if isempty(Y_K3), Y_K3 = 0; end
if isempty(Y_K2), Y_K2 = 0; end
if isempty(Y_K1), Y_K1 = 0; end

X_K1=filter_in;


Y_K1=((Tk^2)*10000*X_K1+20000*(Tk^2)*X_K2+10000*(Tk^2)*X_K3-Y_K2*(20000*(Tk^2)-8)-Y_K3*(10000*(Tk^2)-200*sqrt(2)*Tk+4))/(4+(Tk^2)*10000+200*sqrt(2)*Tk);
X_K3=X_K2;
X_K2=X_K1;
Y_K3=Y_K2;
Y_K2=Y_K1;

filter_out=Y_K1;

五、Matlab仿真验证

在构建了巴特沃兹滤波器后，我们对比三种性能：

一阶巴特沃兹滤波器,二阶巴特沃兹滤波器连续域,二阶巴特沃兹滤波器的离散化形式：

同理仿真采用上述形式方程，截至频率为100Hz,增益为1，首先看时域特性曲线：

可以看出二阶巴特沃兹滤波器连续域,二阶巴特沃兹滤波器的离散化形式重叠，验证了离散化的有效性，其次也看出一阶巴特沃兹滤波器（由于在通带和阻带区间的频率较多），因此响应较快。

其次我们验证一阶和二阶巴特沃兹滤波器的Bode图（幅频，相频特性）：

蓝色为一阶巴特沃兹滤波器，橘色为二阶巴特沃兹滤波器

从上图，可以看出仿真结果和理论分析一致，附代码：

clear
s=tf('s'); 
G1=100/(s+100);
P=bodeoptions;
P.Grid='on';
y=logspace(0,3);
bode(G1,P,y);

s=tf('s'); 
G1=10000/(s^2+100*sqrt(2)*s+10000);
P=bodeoptions;
P.Grid='on';
y=logspace(0,3);
hold on
bode(G1,P,y);

六、结论

优点：相比于另外两种滤波器，巴特沃兹滤波器在passband具有更加flat的特性,在电机控制中的电流滤波环节中可以减少不必要的低频正增益。

缺点：高阶下的巴特沃兹滤波器在实际系统会增大运算量，一般常采用四阶；另外巴特沃兹滤波器的延迟相位非线性，不利于补偿。