数字滤波器的MATLAB设计与DSP上的实现

　　引言

　　随着信息时代和数字世界的到来，数字信号处理已成为今一门极其重要的学科和技术领域。数字信号处理在通信、语音、图像、自动控制、雷达、军事、航空航天、医疗和家用电器等众多领域得到了广泛的应用。在数字信号处理应用中，数字滤波器十分重要并已获得广泛应用。

　　1 数字滤波器的设计

　　1.1 数字滤波器设计的基本步骤

　　数字滤波器根据其冲激响应函数的时域特性，可分为两种，即无限长冲激响应（IIR）滤波器和有限长冲激响应（FIR）滤波器。IIR滤波器的特征是，具有无限持续时间冲激响应。种滤波器一般需要用递归模型来实现，因而有时也称之为递归滤波器。FIR滤波器的冲激响应只能延续一定时间，在工程实际中可以采用递归的方式实现，也可以采用非递归的方式实现。数字滤波器的设计方法有多种，如双线性变换法、窗函数设计法、插值逼近法和Chebyshev逼近法等等。随着MATLAB软件尤其是MATLAB的信号处理工作箱的不断完善，不仅数字滤波器的计算机辅助设计有了可能，而且还可以使设计达到最优化。

　　数字滤波器设计的基本步骤如下： （1）确定指标 在设计一个滤波器之前，必须首先根据工程实际的需要确定滤波器的技术指标。在很多实际应用中，数字滤波器常常被用来实现选频操作。因此，指标的形式一般在频域中给出幅度和相位响应。幅度指标主要以两种方式给出。第一种是绝对指标。它提供对幅度响应函数的要求，一般应用于FIR滤波器的设计。第二种指标是相对指标。它以分贝值的形式给出要求。在工程实际中，这种指标最受欢迎。对于相位响应指标形式，通常希望系统在通频带中人有线性相位。运用线性相位响应指标进行滤波器设计具有如下优点：①只包含实数算法，不涉及复数运算；②不存在延迟失真，只有固定数量的延迟；③长度为N的滤波器（阶数为N-1），计算量为N/2数量级。因此，本文中滤波器的设计就以线性相位FIR滤波器的设计为例。 （2）逼近 确定了技术指标后，就可以建立一个目标的数字滤波器模型。通常采用理想的数字滤波器模型。之后，利用数字滤波器的设计方法，设计出一个实际滤波器模型来逼近给定的目标。 （3）性能分析和计算机仿真 上两步的结果是得到以差分或系统函数或冲激响应描述的滤波器。根据这个描述就可以分析其频率特性和相位特性，以验证设计结果是否满足指标要求；或者利用计算机仿真实现设计的滤波器，再分析滤波结果来判断。

　　1.2 滤波器的MATLAB设计

　　（1）MATLAB MATLAB是一套用于科学计算的可视化高性能语言与软件环境。它集数值分析、矩阵运算、信号处理和图形显示于一体，构成了一个界面友好的用户环境。它的信号处理工具箱包含了各种经典的和现代的数字信号处理技术，是一个非常优秀的算法研究与辅助设计的工具。在设计数字滤波器时，通常采用MATLAB来进行辅助设计和仿真。 （2）FIR滤波器的MATLAB设计 下面以设计线性相位FIR滤波器为例介绍具体的设计方法。 线性相位FIR滤波器通常采用窗函数法设计。窗函数法设计FIR滤波器的基本思想是：根据给定的滤波器技术指标，选择滤波器长度N和窗函数ω（n），使其具有最窄宽度的主瓣和最小的旁瓣。其核心是从给定的频率特性，通过加窗确定有限长单位脉冲响应序列h(n)。工程中常用的窗函数共有6种，即矩形窗、巴特利特（Bartlett）窗、汉宁（Hanning）窗、汉明（Hamming）窗、布莱克曼（Blackman）窗和凯塞（Kaiser）窗。 假设实际工程需要设计一个线性相位带通FIR滤波器指标如下： fn=[1000,1375,3625,4000]; a=[0,1,0]; dev=[0.0005,0.05,0.0005]; 阻带最小衰减60dB，再设采样频率fs=10kHz，则根据阻带最小衰减来选择凯塞（Kaiser）窗，利用MATLAB的求阶函数和FIR滤波器的设计函数，可以快速地设计出所需的数字滤波器。这两个设计函数如下： [N,Wn,beta,ftype]=kaiserord(fn,a,dev,fs)； b=fir1(n,Wn,ftype,Kaiser(n+1,beta)； 最后，利用的滤波器分析函数freqz分析所设计出的滤波器的幅频特性和相频特性，并用图形显示函数plot将它们显示出来，如图1所示。由图1可见，设计结果满足指标要求。

　　2 数字滤波器的实现方法

　　数字滤波器的实现方法一般有以下几种。 ①采用加法器、乘法器、延时器设计专用的滤波电路。 ②在通用计算机系统中加上专用的加速处理机设计实现。 ③用通用的可编程DSP芯片实现。 ④用专用的DSP芯片实现。在一些特殊的场合，要求的信号处理速度极高，用通用DSP芯片很难实现。这种芯片将相应的滤波算法在芯片内部用硬件实现，无需进行编程。 ⑤采用FPGA/CPLD设计实现。 在上述几种方法中，第②种方法的缺点是速度较慢，一般可用于DSP算法的模拟。第①、④种方法专用性强，应用受到很大的限制。第③、⑤种方法都可以通过编程来实现各种数字滤波，但是，第③种因有专用的指令来实现滤波运算编程实现容易，而第⑤种方法编程实现较为困难。

　　3 数字滤波器的DSP实现

　　DSP是一种实时、快速、特别适合于实现各种数字信号处理运算的微处理器。由于它由具有丰富的硬件资源、改进的哈佛结构、高速数据处理能力和强大的指令系统，而在通信、航空、航天、雷达、工业控制、网络及家用电器等各个领域得到广泛应用。DSP分为定点和浮点两种，下面以定点DSP芯片为例，讨论FIR滤波器实现的几个关键问题。所讨论的这些问题，在DSP系统设计中有实际的参考和应用价值。

　　3.1 定点数的定标

　　在定点DSP芯片中，采用定点数进行数值运算，其操作数一般采用整型数来表示。一个整型数的最大表示范围取决于DSP芯片所给定的字长，一般为16位或24位。显然，字长越长，所能表示的数的范围越大，精度也越高。在滤波器的实现过程中，DSP所要处理的数可能是整数，也可能是小数或混合小数；然而，DSP在执行算术运算指令时，并不知道当前所处理的数据是整数还是小数，更不能指出小数点的位置在哪里。因此，在编程时必须指定一个数的小数点处于哪一位，这就是定标。通过定标，可以在16位数的不同位置上确定小数点，从而表示出一个范围大小不同且精度也不同的小数。例如：在Q15中，1080H=0.12890625；在Q0时，1080H=4224。 在使用定点DSP时，如何选择合适的Q值是一个关键性问题。就DSP运算的处理过程来说，实际参与运算的都是变量，有的是未知的，有的则在运算过程中不断改变数值，但它们在一实际工程环境中作为一个物理参量而言都有一定的动态范围。只要个动态范围确定了，Q值也就确定了。因此，在程序设计前，首先要通过细致和严谨的分析，找出参与运算的所有变量的变化范围，充分估计运算中可能出现的各种情况，然后确定采用何种定标标准才能保证运算结果正确可靠。这里，所讨论的理论分析法和统计分析法确定变量绝对值最大值|max|，然后根据|max|再确定Q值。但是，DSP操作过程中的意外情况是无法避免的，即使采用统计分析法也不可能涉及到所有情况。因此，在定点运算过程中应该采取一些判断和保护措施（特别是在定点加法中）。另外，在数字信号处理中的大量运算是乘法和累加，应尽量采用纯整数或纯小数运算，即全部变量都用Q0或Q15格式表示。这样做的好处是操作简单、编程方便。只有当纯整数或纯小数运算不能满足变量的动态范围和精度要求时，才采用混合小数表示法进行定点运算。

　　3.2 误差问题

　　因为在用定点DSP实现时，所有的数据都是定长的，运算也都是定点运算，因而会产生有限字长效应。所产生的误差主要包括：数模转换引起的量化误差、系数量化引起的误差以及运算过程中的舍入误差。在用定点DSP时，产生误差是不能避免的，但是可以通过一些方法减小误差。如，可以用两个存储单元来表示一个数，运算时使用双字运算；可以根据需要要将滤波器系数都用双字表示，也可以只将一半的系数用双字表示，视需要而定。另外，FIR数字滤波器和IIR数字滤波器所引入的量化误差是不一样的。FIR数字滤波器主要采用非递归结构，因而在有限精度的运算中都是稳定的；而IIR数字滤波器是递归结构，极点必须在z平面单位圆内才能稳定，这种结构运算中的四舍五入处理有时会引起寄生振荡。除了有限字长效应以外，不同结构引入的误差也有所不同。在实际设计中，要注意实现中的误差问题。在选择不同的结构时，应考虑它们所引入的误差，并用高级语言进行定点仿真，以比较不同结构下误差的大小，从而作出合理选择。

　　3.3 循环寻址

　　循环寻址（circular addressing）是DSP中经常用到的一种寻址方式。该寻址方法可以对一块特定存储区实现循环的操作。可以把循环寻址理解为实现一个滑动窗，新数据引入后将覆盖老的数据，便得该窗中包含了需处理的最新数据。在数字信号处理中的FIR、卷积等运算中，循环寻址具有极其重要的意义。 在TI的DSP中，循环寻址通过如下方法实现。 ①设定BK（寄存器块大小）值，以确定循环寻址缓冲区的大小，也可将它看作是循环的周期。 ②设定缓冲区的底部地址。必须注意：其低N位为零，其中N为满足式（1）的最小N值。 2 N>BK （1） ③用辅助寄存器间接寻址循环缓冲区。 图2 输入数据频谱图 图3 输出数据频谱图 3.4 实现举例 根据上述设计出来的FIR滤波器，考虑工程实际的需要（精度的要求）和采用定点DSP芯片的方便实现，选择Q15定标。为此，必须对输入数据和滤波器的系数进行归一化处理。输入数据的归一化处理可通过设置A/D转换的参考电平来实现。滤波器的系数归一化只要求取系数中的最大值，再用这个最大值去系数便可。 由于一个N（设N为偶数）阶的FIR滤波器具有系数对称特性，其输出方程可以写为： y(n)=ao[x(n)+x(n-N+1)]+a1[x(n-1)+x(n-N+2)]+…+aN/2-1[X(n-N/2+1)+x(n-N/2)] (2) 根据式（2）可建立如下实现算法： ①在数据存储中开辟两个循环缓冲区，New循环缓冲区存放新数据，O1d循环缓冲区中存放老数据。循环缓冲区的长度为N/2。 ②设置循环缓冲区指针，AR2指向New缓冲区中最新的数据，AR3指向O1d缓冲区中最老的数据。 ③在程序存储器中设置系数表。 ④（AR2）+（AR3）→AH（累加器A的高位）；（AR2）-1→AR2；（AR3）-1→AR3。 ⑤将累加器B清零，重复执行下列操作N/2次：（AH）*系数ai+(B)→B，系数指针（PAR）加1，（AR2）+（AR3）→AH，AR2和AR3减1。 ⑥保存和输出结果（结果在BH中）。 ⑦修正数据指针，让AR2和AR3分别指向New缓冲区最新数据和O1d缓冲区中最老的数据。 ⑧用New缓冲区中最老的数据替代O1d缓冲区中最老的数据。O1d缓冲区指针减1。 ⑨输入一个新数据替代New缓冲区中最老的数据。 重复执行第④～⑨步。 根据上述算法编制程序[4、5]，在CCS5000上进行仿真调试运行，并分析输入数据和输出的频谱，结果如图2、图3所示。由图可见，滤波器实现了目标要求。 最后，将程序移植为闻亭公司的TMS320VC5409评估板的MCBSP0的中断服务程序，并将形成的可执行文件加载到评估板上运行。由TMS320VC5409评估板的模拟输入口输入模拟信号，经TLC320A/D转换芯片，按8kHz的采样频率采样转换成数字信号，输入DSP。滤波后经D/A输出模拟量，结果证实该程序可实现对采样率为8kHz的模拟信号进行实时滤波处理。 结语 数字滤波器的应用十分广泛，运行MATLAB语言，能很容易地设计出具有严格要求（如线性相位等）的滤波器。用定点DSP实现滤波器械要考虑DSP的定标、误差、循环寻址等几个关键问题。

　　文中实例是为了表明，可方便地用DSP实现模拟信号的实时滤波处理，所采用的采样频率并不高。如果DSP采用更高的时钟，它的处理速度将更快，将能够满足更高采样率的数字信号的实时滤波处理。