通用的一阶IIR数字高通滤波器的实现

书接上回

上一篇讲述了一阶IIR滤波器的一些特点，并设计了一个高通的IIR。实际应用中往往需要动态调整滤波器的参数满足不同场合，不同平台的需要，本章基于上一节的高通滤波器来继续给出一个通用的灵活的高通IIR：$$H(s)=\frac{s}{s+{\Omega }_c}$$
数字滤波器的截止频率计算公式：$${\omega }_c=\frac{2\pi f}{f_s}$$

开始推算

数字和模拟频率的转换关系：$${\Omega }_c=\frac{{\omega }_c}{T}$$
直接带入传递函数得到：
$$H(s)=\frac{s}{s+\frac{{\omega }_c}{T}}$$
带入双线性变换公式：$$s=\frac{2}{T}\frac{1-z^{-1}}{1+z^{-1}}和{\omega }_c=\frac{2\pi f}{f_s}$$
我们得到
$$H(s)=\frac{\frac{2}{T}\frac{1-z^{-1}}{1+z^{-1}}}{\frac{2}{T}\frac{1-z^{-1}}{1+z^{-1}}+\frac{\frac{2\pi f}{f_s}}{T}}=\frac{1-z^{-1}}{(\pi \frac{f_c}{f_s}+1)+(\pi \frac{f_c}{f_s}-1)z^{-1}}=\frac{Y(Z)}{X(Z)}$$
根据这个得出时域差分方程$$y(n)=\frac{1}{\pi \frac{f_c}{f_s}+1}(x(n)-x(n-1))-\frac{\pi \frac{f_c}{f_s}-1}{\pi \frac{f_c}{f_s}+1}y(n-1)$$
同理得出预畸变(双线性变换对频率作预畸变处理)：$${\Omega }_c=\frac{2}{T}\tan \frac{{\omega }_c}{2}$$的差分方程$$y(n)=\frac{1}{tan(\frac{\pi f_c}{f_s})+1}(x(n)-x(n-1))-\frac{tan(\frac{\pi f_c}{f_s})-1}{tan(\frac{\pi f_c}{f_s})+1}y(n-1)$$
为何提出这两种形式呢？ 主要是因为在低频段两者基本可以互换，而前者能省掉tan计算，对于一些比较弱鸡的嵌入式系统，三角计算是不支持的。
下边给了一个例子，在48k采样率的情况下，1000Hz以内的截止频率，两者的系数计算结果基本可以互换。

PRE_WARP: temp_alpha1=-0.986995 temp_beta1=0.993497 freqz: 100, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-16170 temp_beta1=16277 freqz: 100, srate:48000
NORMAL: alpha1=-0.986995 beta1=0.993498 freqz: 100, srate:48000
NORMAL TRANSTO INT: a0=-16170 a1=16277 freqz: 100, srate:48000


PRE_WARP: temp_alpha1=-0.974157 temp_beta1=0.987078 freqz: 200, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-15960 temp_beta1=16172 freqz: 200, srate:48000
NORMAL: alpha1=-0.974158 beta1=0.987079 freqz: 200, srate:48000
NORMAL TRANSTO INT: a0=-15960 a1=16172 freqz: 200, srate:48000



PRE_WARP: temp_alpha1=-0.961481 temp_beta1=0.980741 freqz: 300, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-15752 temp_beta1=16068 freqz: 300, srate:48000
NORMAL: alpha1=-0.961486 beta1=0.980743 freqz: 300, srate:48000
NORMAL TRANSTO INT: a0=-15752 a1=16068 freqz: 300, srate:48000


PRE_WARP: temp_alpha1=-0.948965 temp_beta1=0.974482 freqz: 400, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-15547 temp_beta1=15965 freqz: 400, srate:48000
NORMAL: alpha1=-0.948976 beta1=0.974488 freqz: 400, srate:48000
NORMAL TRANSTO INT: a0=-15548 a1=15966 freqz: 400, srate:48000



PRE_WARP: temp_alpha1=-0.936602 temp_beta1=0.968301 freqz: 500, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-15345 temp_beta1=15864 freqz: 500, srate:48000
NORMAL: alpha1=-0.936624 beta1=0.968312 freqz: 500, srate:48000
NORMAL TRANSTO INT: a0=-15345 a1=15864 freqz: 500, srate:48000


PRE_WARP: temp_alpha1=-0.924390 temp_beta1=0.962195 freqz: 600, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-15145 temp_beta1=15764 freqz: 600, srate:48000
NORMAL: alpha1=-0.924428 beta1=0.962214 freqz: 600, srate:48000
NORMAL TRANSTO INT: a0=-15145 a1=15764 freqz: 600, srate:48000


PRE_WARP: temp_alpha1=-0.912326 temp_beta1=0.956163 freqz: 700, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-14947 temp_beta1=15665 freqz: 700, srate:48000
NORMAL: alpha1=-0.912384 beta1=0.956192 freqz: 700, srate:48000
NORMAL TRANSTO INT: a0=-14948 a1=15666 freqz: 700, srate:48000


PRE_WARP: temp_alpha1=-0.900404 temp_beta1=0.950202 freqz: 800, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-14752 temp_beta1=15568 freqz: 800, srate:48000
NORMAL: alpha1=-0.900491 beta1=0.950245 freqz: 800, srate:48000
NORMAL TRANSTO INT: a0=-14753 a1=15568 freqz: 800, srate:48000


PRE_WARP: temp_alpha1=-0.888622 temp_beta1=0.944311 freqz: 900, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-14559 temp_beta1=15471 freqz: 900, srate:48000
NORMAL: alpha1=-0.888744 beta1=0.944372 freqz: 900, srate:48000
NORMAL TRANSTO INT: a0=-14561 a1=15472 freqz: 900, srate:48000


PRE_WARP: temp_alpha1=-0.876976 temp_beta1=0.938488 freqz: 1000, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-14368 temp_beta1=15376 freqz: 1000, srate:48000
NORMAL: alpha1=-0.877141 beta1=0.938571 freqz: 1000, srate:48000
NORMAL TRANSTO INT: a0=-14371 a1=15377 freqz: 1000, srate:48000



PRE_WARP: temp_alpha1=-0.865464 temp_beta1=0.932732 freqz: 1100, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-14179 temp_beta1=15281 freqz: 1100, srate:48000
NORMAL: alpha1=-0.865681 beta1=0.932840 freqz: 1100, srate:48000
NORMAL TRANSTO INT: a0=-14183 a1=15283 freqz: 1100, srate:48000



PRE_WARP: temp_alpha1=-0.854081 temp_beta1=0.927040 freqz: 1200, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-13993 temp_beta1=15188 freqz: 1200, srate:48000
NORMAL: alpha1=-0.854359 beta1=0.927179 freqz: 1200, srate:48000
NORMAL TRANSTO INT: a0=-13997 a1=15190 freqz: 1200, srate:48000


PRE_WARP: temp_alpha1=-0.842824 temp_beta1=0.921412 freqz: 1300, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-13808 temp_beta1=15096 freqz: 1300, srate:48000
NORMAL: alpha1=-0.843174 beta1=0.921587 freqz: 1300, srate:48000
NORMAL TRANSTO INT: a0=-13814 a1=15099 freqz: 1300, srate:48000



PRE_WARP: temp_alpha1=-0.831691 temp_beta1=0.915846 freqz: 1400, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-13626 temp_beta1=15005 freqz: 1400, srate:48000
NORMAL: alpha1=-0.832123 beta1=0.916061 freqz: 1400, srate:48000
NORMAL TRANSTO INT: a0=-13633 a1=15008 freqz: 1400, srate:48000


PRE_WARP: temp_alpha1=-0.820679 temp_beta1=0.910339 freqz: 1500, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-13446 temp_beta1=14915 freqz: 1500, srate:48000
NORMAL: alpha1=-0.821204 beta1=0.910602 freqz: 1500, srate:48000
NORMAL TRANSTO INT: a0=-13454 a1=14919 freqz: 1500, srate:48000



PRE_WARP: temp_alpha1=-0.809784 temp_beta1=0.904892 freqz: 1600, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-13267 temp_beta1=14825 freqz: 1600, srate:48000
NORMAL: alpha1=-0.810414 beta1=0.905207 freqz: 1600, srate:48000
NORMAL TRANSTO INT: a0=-13277 a1=14830 freqz: 1600, srate:48000


PRE_WARP: temp_alpha1=-0.799004 temp_beta1=0.899502 freqz: 1700, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-13090 temp_beta1=14737 freqz: 1700, srate:48000
NORMAL: alpha1=-0.799751 beta1=0.899876 freqz: 1700, srate:48000
NORMAL TRANSTO INT: a0=-13103 a1=14743 freqz: 1700, srate:48000


PRE_WARP: temp_alpha1=-0.788336 temp_beta1=0.894168 freqz: 1800, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-12916 temp_beta1=14650 freqz: 1800, srate:48000
NORMAL: alpha1=-0.789213 beta1=0.894607 freqz: 1800, srate:48000
NORMAL TRANSTO INT: a0=-12930 a1=14657 freqz: 1800, srate:48000



PRE_WARP: temp_alpha1=-0.777778 temp_beta1=0.888889 freqz: 1900, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-12743 temp_beta1=14563 freqz: 1900, srate:48000
NORMAL: alpha1=-0.778798 beta1=0.889399 freqz: 1900, srate:48000
NORMAL TRANSTO INT: a0=-12759 a1=14571 freqz: 1900, srate:48000


PRE_WARP: temp_alpha1=-0.767327 temp_beta1=0.883663 freqz: 2000, srate:48000
PRE_WARP: TRANSTO INT: temp_alpha1=-12571 temp_beta1=14477 freqz: 2000, srate:48000
NORMAL: alpha1=-0.768503 beta1=0.884252 freqz: 2000, srate:48000
NORMAL TRANSTO INT: a0=-12591 a1=14487 freqz: 2000, srate:48000

浮点转定点

低端的数字信号处理器大都是定点处理器，所以真要这个滤波器活蹦乱掉的跑起来，还得做一个浮点转定点的工作。听着很高大上，其实很简单，就是系数转换为Qxx的整形格式（以2的整数次幂放大），然后一顿加乘计算后，再等量缩减。细心的读者可以发现脚本例子中“TRANSTO INT”一行的系数算出了一个巨大的整数值，这个值就是浮点转定点的系数。实在觉得费尽，可以找一些现成的代码，下面的代码就是一个Q15的转换函数。

int To_Q15_Fixed_Point(double Data_Float)
{
    unsigned int  Data_Q15_Format = 0;

    if(Data_Float == 1) return (unsigned short)0x7FFF;
    else if(Data_Float == -1) return (unsigned short)0x8000;
    else if(Data_Float < 0)
    {
        Data_Float = -Data_Float;
        Data_Q15_Format |= (unsigned short)0x8000;
    }

    Data_Q15_Format |= (short)(Data_Float*0x8000);

    if(Data_Q15_Format&0x8000)
    {
        Data_Q15_Format &= 0x7FFF;
        Data_Q15_Format = ((~Data_Q15_Format) + 0x0001);
    }

    return (int)Data_Q15_Format;
}

最后

有了系数和差分方程，写函数这件事就是最简单的了。在上一篇的基础上，我们彻底解放了这个一阶iir的高通滤波器，使用的时候需要评估采用哪种系数计算方法更满足要求。