经典基于dsp数字信号处理器的fir滤波器的设计和实现

经典基于dsp数字信号处理器的fir滤波器的设计和实现内容摘要：

收端都是发出或接收一路信号。 所以，前两种滤波器还是现在滤波器设计的主要方面，例如在线谱分析、基音检测、线性预测编码等方面都有着广泛的应用。 本来在计算量相等的情况下， IIR 数字滤波器比 FIR 滤波器的幅频特性优越，频率选择性也好，但是，它有着致命的缺点，相位特性不好控制。 它的相位特性 f(w) argH( jwe )是使频率产生严重的非线性的原因，这 种‘与’的非线性关系，使数字滤波器与模拟滤波器在响应与频率的对应关系上发生了畸变。 如果需要线性相位，就必须用全通网络进行复杂的相位校正，但是，在对程序运行周期数要求十分严格的 DSP 处理中加上一个全通均衡器是十分浪费资源的。 另外，即使加上全通均衡器，对于因果的 IIR 滤波器，仍将得不到线性的相位。 在现代电子系统中，如图像处理、数据传输等波形传递系统中都越来越多的要求信道具有线性的相位特性。 在这方面 FIR滤波器具有独到的优点，它可以在幅度特性随意设计的同时，保证精确、严格的线性相位，因此这类滤波器应用很广泛。 FIR 滤波器的基本结构 数字滤波是将输入的信号序列，按规定的算法进行处理，从而得到所期望的输出序列。 一个线性位移不变系统的输出序列 y(n)和输入 x(n)之间的关系，应满足常系数线性差分方程，见公式 ，         Mi iNi i inyainxbny 110 0n () 其中， x(n)为输入序列， y(n)为输出序列， kk ba和 为滤波器系数， N是滤波器的阶数。 若上式中所有的 kb 均为零，则有 FIR滤波器的差分方程为：   10 )()( Nk k knxany () 对上式进行 Z变换得到 FIR滤波器的传递函数为：        10Ni kk zbzX zYzH () 由上式可以看出， H(z)是 1Z 的 N1次多项式，它在 z平面内有 N1个零点，同时在原点处有 N1个重极点。 N阶滤波器通常采用 N个延迟单元、 N个加法器与 N+1个乘法器 ，取图 31 中 (a)、 (b)两种结构。 图 FIR滤波器的一般结构 因为 FIR滤波器的单位抽样响应是有限长的，所以它永远是稳定的。 另外，若对 h(n)提出一些约束条件，那么可以很容易地使 H(z)具有线性相位，这在信号处理的很多领域是非常重要的。 FIR滤波器的设计任务，是要决定一个转移函数 H(z)，使它的频率响应满足给定的要求。 这里所说的要求，除了通带频率 p 、阻带频率及两个带上的最大和最小衰减 p 和 s 外，很重要的一条是保证 H(z)具有线性相位。 FIR 滤波器的常规设计方法 FIR滤波器的设计任务就是给定要求的频率特性，按一定的最佳逼近准则，选取滤波器转移函数 H(z)中的各个参数 h(n)，即滤波器的单位抽样响应及阶数 N，使得频率特性满足设计要求。 通常 FIR滤波器的设计方法主要有三种：窗函数法、频率抽样法和切比雪夫等波纹逼近法。 其中窗函数法可以应用比较现成的窗函数，因而设计简单，在指标要求不高的场合使用方便灵活。 下面我们来简单介绍一下这三种设计方法。 窗函数法 窗函数法也称为傅立叶级数法。 理想的数字滤波器频率特性 )( jweH 是无法实现的， FIR的设计就是要寻找一个可以得到的频率特性 )( jweH =10 )(Nnjwnenh 来逼近 )( jweH ，这相当于用一个可实现的单位脉冲响应 h(n)去逼近一个理想单位脉冲响应 )(nhd。 )(nhd 可由理想频率特性 )( jwd eH 通过傅氏反变换得到，   deeHnh jwjwdd  )(21)( () 一般来说，这样得到的理想单位脉冲响应序列 )(nhd 是个无限长序列，因而是非因果的。 设有一个截止频率为 c 的理想线性相位低通，延时为 τ，其频率特性是：       ccjw ajwd eeH 0 0 () 得到： )(nhd ＝  )( )(sin   n nc  n () 这是一个以ｎ＝τ为中心偶对称的无限长非因果序列，要想用一个有限长的因果序列去逼近它，最简单的方法是截取 ｎ 从 0到 N1的一段来表示它，即 h(n)= )(nhd )10(  Nn ；其他 N： h(n)=0。 同时，为了保证线性相位，还要满足偶对称 h(n)=h(N1n)。 这就好像 通过一个窗口观看到的一段 )(nhd ，因此 h(n)就表示成 )(nhd 和一个窗口函数的乘积，这样对 h(n)的求解就变为 h(n)＝ )(nhd * nW ，这里的 nW 就称为窗口函数 , 既然一个频域上的标准的矩形窗口对应于时域是一个无限长的序列 , 那么在时域上截取一段势必造成频域的矩形窗口的失真。 结果就是截取出的信号也相应失真，为了补 偿这种失真，只有改变原来窗口的形状，修正经过时域截取后的窗口失真。 窗函数设计方法的基本步骤是： (1) 把  jwd eH 展成 FS，得 )(nhd ； (2) 对 )(nhd 自然截短到所需的长度，如 2M+1； (3) 将截短后的 )(nhd 右移 M 个采样间隔，得 h(n)； (4) 将 h(n)乘以合适的窗口，即得所要滤波器的冲击响应，窗函数以 n=M 对称。 利用所求得的单位抽样响应，即 可用硬件构成滤波器的转移函数 H(z)，也可利用 h(n)在计算机上用软件来实现滤波。 频率抽样法 窗函数法是从时域出发，用窗函数截取理想的 nhd 得到 h(n)，以此有限长的 h(n)近似 nhd ，这样得到的频率响应  jeH 逼近于理想的频响  jwd eH。 频率抽样法是从频率出发， 将给定的理想频响  jwd eH 加以等间隔抽样。  jwd eH  kHdNk  2 () 然后以此 kHd 作为 FIR滤波器的频率响应抽样值 H(k)，再根据 DFT（离散付氏变换）定义由频域这 N个抽样值来唯一确定一个有限长序列 h(n)，同样也可以算出 FIR滤波器的系统函数 H(z)及频率响应  jeH ，可以推出频率响应  jeH 是频率抽样值 H(k)与线性相位因子   2/1 Nje  及如下内插函数 S(ω , k) 的线性组合。 S(ω , k)= NkNeN Nkj2sin2sin1 () 所以，在各频率取样点上，实际滤波器的频响是严格地和所要求的滤波器的频响一致的，逼近误差为零，但在抽样点之间的频响是各取样点的内插函数的延伸叠加而成，有一定的逼近误差，误差大小取决于频率响应曲线的圆滑程度和抽样点的密度 为了减少误差 就要增加抽样点数目即增大采样频率，抽样点之间的理想频率特性变化越陡，则逼近误 差越大，在理想频率特性的不连续点附近会产生肩峰和纹波。 频率抽样法的优点是可以直接在频域设计，适于利用最优化方法，而且这种方法特别适用于窄带选频滤波器， 但频率抽样法的抽样频率只能是 2π / N 的整数倍或 2π / N 的整数倍加上π / N不能保证截止频率ω c的准确取值，要实现精确的ω c就必须取 N大，相应的计算量也大。 此外，它的阻带最大衰减一般，也只有 3050dB左右， 很难满足频域特性要求较高的场合。 Chebyshev 逼近法 窗函数法和频率采样法设计出的滤波器的频率特性都是在不同意义上对所给 理想频率特性  jwd eH 的逼近。 由数值逼近理论可知，对某个函数 f(x)的逼近一般有以下三种方法： 插值法 (Interpolating Way) 最小平方逼近法 (Least Square Approaching Way) 一致逼近法 (Consistent Approaching Way) 切比雪夫最佳一致逼近的基本思想是，对于给定区间 [a， b]上的连续函数 xf ，在所有 n次多项式的集合 n 中，寻 找一个多项式 p (x)，使它在 [a， b]上对 xf 的偏差和其它一切属于 n 的多项式 p(x)对 f(x)的偏差相比是最小的，即          xfxpxfxp  m a xm i nˆm a x () 切比雪夫逼近理论，这样的多项式是存在的，且是唯一的，并指出了构造这种最佳一致逼近多项式的方法，就是有名的“交错点组定理”。 切比雪夫逼近理论解决了 p(x)的存在性、唯一性和如何构造等问题。 、 、 等人应用切比雪夫逼近理论提出了一种设计FIR滤波器的计算机辅助算法。 这种算法由于是在一致意义上对  jwd eH 作最佳逼近，因而获得了较好的通带和阻带性能，并能准确地指定通带和阻带的边缘。 但它的效率依赖于初始极值频率点的估计，且通带和阻带内波纹数较多，这是 Chebyshev方法的两个主要缺点。 FIR 滤波器的 MATLAB 实现 MATLAB信号处理工具箱提供了基于窗函数法的 FIR滤波器的设计函数 fir1。 fir1是采用经典窗函数法 设计线性相位 FIR数字滤波器，且具有标准低通、带通、高通和带阻等类型 ]12~11[。 语法格式： B=fir1（ n, nW ） B=fir1（ n, nW ,‘ftype‘） B=fir1（ n, nW ,window） B=fir1（ n, nW ,‘ftype‘， window） 其中， n 为 FIR 滤波器的阶数，对于高通、带阻 滤波器 n 取偶数。 nW 为滤波器截止频率，取值范围为 0~1。 对于带通、带阻滤波器， nW =[ 1W ， 2W ]，且 1W 2W。 ‘ftype‘为滤波器类型。 缺省时为低通或带通滤波器，为‘ high’时是高通滤波器，为‘ stop’时是带阻滤波器。 Window 为窗函数，列向量，其 长度为 n+1；缺省时，自动取 hamming窗。 输出参数 B 为 FIR 滤波器系数向量，长度为 n+1。 带通滤波器的 MATLAB 实现 使用矩形窗、 Hanning窗、 Hamming窗、布莱克曼窗四种窗对带通原型进行截取。 带通滤波器的指标性能给出如下： 下阻带边缘：  s ， dBAs 60 下通带边缘：  p ， dBRp 1 上通带边缘：  p ， dBRp 1 上阻带边缘：  s ， dBAs 60 设计结果如 、 、 、。 图 矩形窗设计 图 图 hamming 窗设计 图 布莱克曼窗设计 低通滤波器的 MATLAB 实现 用窗函数法设计 FIR低通滤波器，其技术指标如下： 250fs kHz, 20pf kHz， 30sf kHz， 通带最大衰减 3dBpA  ，阻带最小衰减50dBsA 。 将其换算成数字域的性能指标如下： 通带截止频率 p，通带最大衰减 3dBpA ； 阻带截止频率  ，阻带最小衰减 50dBsA 。 根据窗函数法的设计原则，由表 ，海明窗（ hamming） 可提供大于 50dB的衰减。 要求滤波器的过渡带为： 0. 24 0. 16 0. 08sp           由表 ，利用海明窗设计的滤波器的过渡带 8/N ，所以低通滤波器单位脉冲响应的长度为： 88 1000 .0 8N    ，取 N=101。 通过 FIR1函数设计出滤波器的单位冲激响应序列，用它来作为在下一章中 DSP设计程序中的系 数。 其设计结果如图。 (a)幅度响应 (b)相位响应 图 hamming窗滤波器的幅度响应与相位响应 4 FIR 滤波器的应用及其 DSP 实现 FFT/IFFT 算法程序及应用 FFT 设计方法 FFT是 DFT的一个快速算法，是为减少 DFT计算次数的一种快速有效的算法。 它是将 DFT分解开来进行运算，理论上是一致的，只是通过分解 DFT运算来达到减少运算量的目的。 其突出的优点在于能快速高效地和比较精确地完成 DFT的计算。 利用一定的运算结构变换，将 N点的 DFT 转化成多个小的点数 DFT 的运算，再利用 knNW 的周期性和对称性，就能大大减少计算量。 FFT 算法将长序列的 DFT分解为短序列 DFT， N 点的DFT先分解为 2个 N/2点的 DFT，每个 N/2点的 DFT又分解为 N/4点的 DFT，如此这般。 这里最小的变换点数即所谓的基数（ radix），因此，基数为 2的 FFT算法的最小变换或称蝶形变换就是 2点 DFT，是最基本的运算单位，一般 N点 FFT对应于 N个输入样值，有 N个频域样值与之对应。 DFT 分解法基本上分为两类：一类是将时间序列 x(。