用Matlab及C语言实现低通滤波器的设计

已知：

低通滤波器的截止频率 fl = 35Hz；即：通带边缘为35Hz，设定通带纹波1db；

采样频率 fs = 400Hz；

设定阻带边缘为 44Hz，衰减为40DB；

Matlab编程：

% 35Hz低通滤波器

Fsam = 400; %采样频率

fp = 31; %通带边缘

Rp = 1; %通带纹波

fs = 44; %阻带边缘

As = 15; %阻带衰减

wp = 2*fp/Fsam; %对奈奎斯特频率（fsam/2）归一化

ws = 2*fs/Fsam; %对奈奎斯特频率（fsam/2）归一化

[n,Wn]= buttord(wp,ws,Rp,As);

[b,a] = butter(n,Wn)

n %输出滤波器的阶数

Wn = Wn*Fsam/2 %将归一化的截止频率还原，单位为Hz；

[X,w] = freqz(b,a,512,Fsam); %求取系统频率响应

plot(w,abs(X)); %画解卷绕后的幅频响应

title('Butterworth Lowpass Filter');ylabel('幅度');xlabel('频率(Hz)');

这个滤波器的阶数为n = 7； 实际截止频率 = Wn*Fsam/2
= 34.9915Hz ；

注意：用matlab运算结果得到的Wn，这个值就像wp，ws一样，是归一化后的值；实际值要乘以奈奎斯特频率Fsam/2（200Hz）；

b =

0.0000 0.0003 0.0008 0.0014 0.0014 0.0008 0.0003 0.0000

a =

1.0000 -4.5340 9.1151 -10.4545 7.3548 -3.1636 0.7685 -0.0812

n = 7

Wn = 34.9915

C语言实现：

/************************* IIR Lowpass Filter **********************************/

//N = 8; -3.13db@35MHz;-22.5db@50Hz;
Rp<1db;

float b[8] = {0.0000,0.0003,0.0008,0.0014,0.0014,0.0008,0.0003,0.0000};

float a[8] = {1.0000,-4.5340,9.1151,-10.4545,7.3548,-3.1636,0.7685,-0.0812};

/*******************************************************************************/

float xBuf1[8];

float yBuf1[8];

float IIR35HzLP(float x)

{

int i;

//运算之前Buf向前移动一个位置，以保存之前Buf的数据；

for(i=7; i>0; i--)

{

yBuf1[i] = yBuf1[i-1]; xBuf1[i] = xBuf1[i-1];

}

xBuf1[0] = x;

yBuf1[0] = 0;

for(i=1;i<8;i++)

{

yBuf1[0] = yBuf1[0] + b[i]*xBuf1[i];//这里有相同系数，可合并来提前计算效率；

yBuf1[0] = yBuf1[0] - a[i]*yBuf1[i];

}

yBuf1[0] = yBuf1[0] + b[0]*xBuf1[0];

return yBuf1[0];

}

函数说明：

进一个x，出一个y；

x对应差分方程的x(n),在函数中的位置为xBuf[0];n为当前最大下标；

x(n-1)对应xBuf[1]...

yBuf[0]对应y(n),

y(n-1)对应yBuf[1]...