软件滤波算法的频率响应

本文从一种软件滤波算法切入,深入理解了滤波算法的频率响应。

引子

单片机大多资源小,算法占用的资源越小越好,有一种占用很小资源的算法,具体公式如下:
$$ y[n] =y[n-1]- {y[n-1] \over m}+x[n]$$
其中:x[n]为采样值,y[n]为滤波后的值,y[n-1]为上一次滤波值,m是放电比例(最好选2的幂次数,可用单片机移位)初始化时如果是跟踪一段时间后使用,可以是任何值,否则可以用采样值乘m初始化。

推导

如果写成这样
$$ y[n] = x[n] + (1- {1 \over m}) * y[n-1]$$ (1)
然后令 $y[n]=m * u[n]$,得
$$ n * u[n] = x[n] + (1- {1 \over m}) * m * u[n-1]$$
两边同时除m,得
$$ u[n] = {1 \over m} * x[n] + (1- {1 \over m}) * u[n-1]$$
令 $ \alpha = {1 \over m}$,得
$$ u[n] = \alpha * x[n] + (1-\alpha) * u[n-1]$$
这是典型的无限脉冲滤波器(IIR),也是最简单的低通滤波器。

如果直接由式(1)推导z变换
$$ Y(z) = X(z) + (1 - {1 \over m}) * z^{-1}Y(z)$$
$$ H(z) = {Y(z) \over {X(z)}} = {1 \over {1-(1-{1 \over m}})z^{-1}}$$
令$a={1-{1\over m}}$,
$$ H(z)={1\over{1-az^{-1}} }$$

IIR vs FIR

IIR(Infinit Impulse Response)滤波器的单位脉冲响应是无限长的,如例子中$H(z)={1\over{1-az^{-1}} }$时间域的序列是
$$ y[n] = x[n] + ay[n-1]$$
根据单位脉冲响应的定义,即$x[n]= \delta[n]=\begin{cases} {1, n=0}\\{0, n \ne 0}\end{cases}$时的y[n]:
$$ y[0] = x[0] + ay[-1] = 1, if n<0 let y[0] = 0$$
$$ y[1] = x[1] + ay[0] = 0 + a = a $$
$$ y[2] = x[2] + ay[1] = 0 + a * a = a^2 $$

$$ y[n] = x[n] + ay[n-1] = a^n $$
所以其单位脉冲响应为$h[n]=a^n u[n]$,其中$u[n]=\begin{cases} {1, n \ge 0 }\\{0, n<0}\end{cases}$,可以看到$n\ge 0$时,h[n]每一项都不为0,也就是无限脉冲响应的含义.IIR滤波器时间域的输出,由输入和之前的输出共同决定,是有反馈的。

FIR(Finit Impulse Response)滤波器的单位脉冲响应是有限长的,其时间域的输出只有输入决定,是无反馈的。如,滑动平均滤波器,时间域$y[n]={1 \over 3}x[n]+{1 \over 3}x[n-1]+{1 \over 3}x[n-2]$,z变换为$H(z)={1\over3}+{1\over3}z^{-1}+{1\over3}z^{-2}$,单位脉冲响应为$h[n]={1\over3}\delta[n]+{1\over3}\delta[n-1]+{1\over3}\delta[n-2]$,其中$\delta[n]=\begin{cases} {1, n=0}\\{0, n \ne 0}\end{cases}$。

频率响应

时间常数

如果将$y[n]=\alpha x[n]+(1-\alpha)y[n-1]$视作RC低通滤波器
250px-1st_Order_Lowpass_Filter_RC.svg.png
$$v_{in}(t)-v_{out}(t)=RC{dv_{out} \over dt}$$

设$v_{in}$序列为($x_{1}, x_{2}, …, x_{n}$), $v_{out}$序列为($y_{1}, y_{2}, …, y_{n}$).上面的公式近似于
$$x_{i} - y_{i} = RC {y_{i}-y_{i-1}}\over{\Delta T}$$

$$ y_{i} = x_{i}({\Delta T}\over {RC+\Delta T}) + y_{i-1}({RC}\over{RC+\Delta T})$$

令 $\alpha = \frac{\Delta T}{RC+\Delta T}$,
$$ y_{i} = {\alpha x_{i} + (1-\alpha) y_{i-1}}$$
其中, $\alpha = \frac{\Delta T}{RC+\Delta T}$
可以得出时间常数$RC={\Delta T}\frac{1-\alpha} {\alpha}$.
考虑到$f_c= {1 \over{2\pi RC}}$,
$$\alpha = \frac{2 \pi \Delta T f_c}{2\pi \Delta T f_c+1}$$
$$f_c=\frac{\alpha} {(1-\alpha)2\pi\Delta T}$$

Reference:

  1. 一个非常适合单片机的滤波算法
  2. Low pass filter in wiki
  3. 指数平滑
  4. IIR