关于MATLAB辅助设计数字滤波器

关于matlab如何辅助设计滤波器最近大致搞明白一些，这里主要是一些比较传统的数字滤波器设计方法（还有一些直接转换到CCS，VHDL或者C Header的这里没有提到），通过matlab算出相应的传递函数，然后C语言进行仿真测试，最后再改写到所需要用的设备上。

首先谈谈matlab的设计方式，主要两种

1.辅助设计 Filter Designe & analysis Tool,通过辅助工具设计,比较省事

进去观察下有哪些选项可以设置

大部分还是基本看了就知道需要怎么设，然后就是准备导出传递函数

FILE中的EXPORT,然后设置10进制输出

最后就能导出相应的传函

%% Generated by MATLAB(R) 7.8 and the Signal Processing Toolbox 6.11.%% Generated on: 06-Mar-2012 23:13:19%% Coefficient Format: Decimal% Discrete-Time IIR Filter (real)                                     % -------------------------------                                     % Filter Structure    : Direct-Form II, Second-Order Sections         % Number of Sections  : 2                                             % Stable              : Yes                                           % Linear Phase        : No                                                                                                                 SOS matrix:                                                           1  -0.3274034970481377  1  1  -1.487740666582311   0.623662226459985921  -1.4109006645592876  1  1  -1.5591910986997555  0.90887414412475032                                                                      Scale Values:                                                         0.072426444477067295                                                  0.59358927160118591   

得出一个SOS矩阵

这个可以通过下图看下怎么转换

对于以零点-极点-增益为参数的传递函数H（Z）,可用函数zp2sos将其转换为LX6sos数组。函数sos2ss将二阶滤波器节映射为直接II型状态空间形式，而函数sos2tf将一组  2阶滤波器节转换为传递函数。对于滤波器节为1阶的情况，系数b2i和a2i将被没为零。函数sos2zp将2阶滤波器节转换为零点－极点－增益的形式。

2.通过matlab函数进行实现

由于设计的方法挺多了，就不列举了，直接通过函数调用就行了

这里主要是为了方便实现转换为比较容易实现的是级联的设计方式

相关的函数

function [b0,B,A] = dir2cas(b,a);% 直接型到级联型的型式转换(复数对型)% ---------------------------------------------------------% [b0,B,A] = dir2cas(b,a)%  b = 直接型的分子多项式系数%  a = 直接型的分母多项式系数% b0 = 增益系数%  B = 包含各bk的K乘3维实系数矩阵%  A = 包含各ak的K乘3维实系数矩阵% compute gain coefficient b0b0 = b(1); b = b/b0;a0 = a(1); a = a/a0;b0 = b0/a0;%M = length(b); N = length(a);if N > Mb = [b zeros(1,N-M)];elseif M > Na = [a zeros(1,M-N)]; N = M;elseNM = 0;end%K = floor(N/2); B = zeros(K,3); A = zeros(K,3);if K*2 == N;b = [b 0];a = [a 0];end%        broots = cplxpair(roots(b));aroots = cplxpair(roots(a));for i=1:2:2*KBrow = broots(i:1:i+1,:);Brow = real(poly(Brow));B(fix((i+1)/2),:) = Brow;Arow = aroots(i:1:i+1,:);Arow = real(poly(Arow));A(fix((i+1)/2),:) = Arow;end     获得了相关的传函之后，就可以进行测试了，测试也是两部分首先通过matlab测试一些相关的幅频特性%对于有理分式H(jw),MATLAB提供freqs函数处理方法。调用格式为H=freqs(b,a,w)%其中，b为分子多项式的系数，a为分母多项式系数，w为需计算的频率特性clear allclcw=linspace(0,5,200);b=[1];a=[1 2 2 1];H=freqs(b,a,w);plot(w,abs(H));

测试下

其次是用C语言大致模拟下，大致也体现了怎么使用传递函数，如下

#include <math.h>#include <stdio.h>#include <string.h>#include <stdlib.h>#include <conio.h>#include <dos.h>#include <time.h>#include <io.h>#include <process.h>void main(){  int i;  int x0,x1,x2,y0,y1,y2,q0,q1,q2;  FILE *fp;  x0=0;x1=0;x2=0;  y0=0;y1=0;y2=0;  q0=0;q1=0;q2=0;  system("cls");  fp=fopen("C:\\data0.txt","wt");  printf("Press any key to continue when ready\n");  printf("or Press ESE to Cancel\n");  getch();  for(i=0;i<500;i++)  {    x2=1600;    if(i==200) x2=3200;    if(i==350) x2=800;    y2=(113*x2-107*x1+113*x0+1618*y1-662*y0)/1024;    q2=(252*y2-436*y1+252*y0+1818*q1-913*q0)/1024;    y0=y1;y1=y2;x0=x1;x1=x2;q0=q1;q1=q2;    printf("%d %d\n",x2/16,q2/16);    fprintf(fp,"%d %d\n",x2/16,q2/16);  }  fclose(fp);  system("pause");}

这里113,107,113,1618,662为传函的参数调理之后，如下

这样大致就设计好了，这个C语言的版本只是一个直流信号的测试，不过大致可以明白，将传函的微分方程改写成差分方程，并进行计算，通过这个思想，可以很方便的将其移植到其他的单片机或者FPGA上，还是比较好实现的。