基于FPGA的数字滤波器设计（IIR滤波）

 基本原理

1. IIR数字滤波器设计的基本原理

基本原理和结构。IIR滤波器，即无线脉冲响应滤波器，其答案为脉冲响应是无限长的，传递函数可以表示为式。

IIR滤波器有直接I型，直接II型，级联型及并 联型4种常用的结构形式，其中级联型结构便于实现，且受参数量化影响较小，因此使用较为广泛。由差分方程可得，输出信号由两部分组成：第 一部分∑M i＝0 x(n-i)b(i)表示将输入信号进行延时，组成M节延时网络，相当于FIR滤波器得横向网络，实现系统的零点。第二部分∑N l＝1y(n-l)a(l)表示将输出信号进行延时，组成N节点的延时网络，每节延时抽头后与常数相乘，并将乘法结果相加。由于这部分是对输出的延时，故为反馈网络，实现系统的极点。信号流图如图１所示。

图 1 直接Ｉ型IIR滤波器信号流图

对于N阶差分方程，直接II型结构只需要N个延 时单元，比直接I型结构延时单元少一半，所以在软件实现时可接收寄存器，相比直接I型结构具有明显的优势。

2设计内容

本设计中有三种正弦信号混合而成，分别是5*sin(2*pi*200*t)， 5*sin(2*pi*1000*t)， 5*sin(2*pi*2000*t)，要保留的是5*sin(2*pi*200*t)正弦信号。然后将三种信号单独的时域图和混合信号的时域图展示出来，然后再通过IIR滤波，可得到滤波之后的三个信号，可以看出和原信号基本吻合，以下是matlab实现源码：

第一步，在matlab中生成x = 5*sin(2*pi*200*t) + 5*sin(2*pi*1000*t) + 5*sin(2*pi*2000*t)的混频信号，其程序如下：

clear;%清屏
clc;
fs = 8000;N=1024;n=0:N-1;t=n/fs;%sin信号的参数
x = 5*sin(2*pi*200*t) + 5*sin(2*pi*1000*t) + 5*sin(2*pi*2000*t);%混频信号
% c = min(x); %得出信号的最小值，以便知道怎么取正
f = round(x+14)';%根据所得最小值，把信号取正，以便FPGA采用
%c = max(f); %可以得出取正之后数据的最大值，从而知道ROM当中的数据长度
fid = fopen('b.mif','wt');  %将所得的数据存在b.mif当中，（数据文件类型也可以是txt类型）
fprintf(fid,'%g\n',f);
fclose(fid);
[a,b] = size(f); %可以得出数据的位数，好设置ROM的容量
I = a*b;

第二步，产生原始信号各分量的时域图及原始混频信号的时域频域图，然后对原始信号进行低通高通带通带阻的滤波，最后将滤波后的信号的时域频域输出。

clear;%清屏
clc;
fs = 8000;N=1024;n=0:N-1;t=n/fs;%sin信号的参数
x = 5*sin(2*pi*200*t) + 5*sin(2*pi*1000*t) + 5*sin(2*pi*2000*t);%混频信号
figure(1);%原始混频信号
subplot(4,1,1);plot(t,sin(2*pi*200*t));
title('x=sin(2*pi*200*t)');
subplot(4,1,2);plot(t,sin(2*pi*1000*t));
title('x=sin(2*pi*1000*t)');
subplot(4,1,3);plot(t,sin(2*pi*2000*t));
title('x=sin(2*pi*2000*t)');
subplot(4,1,4);plot(t,x);
title('原始混频信号（时域）');
%频谱分析
y = (fft(x,N));%做fft变换
Y2 = (y/N);
Y1 = Y2(1:N/2+1);%此时选取前半部分，fft之后为对称的双边谱
Y1(2:end-1) = 2*Y1(2:end-1);
f = fs*(0:(N/2))/N;
figure(2);
plot(f,abs(Y1));
title('原始混频信号（频谱）');
%低通滤波
Wp = 200/4000;%0.05
Ws = 400/4000;%0.1
[n,Wn] = buttord(Wp,Ws,1,20); %3,60
[b,a] = butter(n,Wn,'low');
 c = min(b);   %将滤波参数取正，以便FPGA当中使用
B = round(b*20000)';
A = round(a*20000)';
c = max(A);
fid = fopen('B.txt','wt');  %将所得的数据存在b.mif当中，（数据文件类型也可以是txt类型）
fprintf(fid,'%g\n',B);
fclose(fid);
fid = fopen('A.txt','wt');  %将所得的数据存在b.mif当中，（数据文件类型也可以是txt类型）
fprintf(fid,'%g\n',A);
fclose(fid);
z = filter(b,a,x);
Z = fft(z,N);
Z2 = (Z/N);
Z1 = Z2(1:N/2+1);%此时选取前半部分，fft之后为对称的双边谱
Z1(2:end-1) = 2*Z1(2:end-1);
figure(3);
subplot(2,1,1);plot(t,z);
title('低通滤波信号（时域）');
subplot(2,1,2);plot(f,Z1);
title('低通滤波信号（频域）');
%高通滤波
Wp = 2000/4000;
Ws = 1000/4000;[n,Wn] = buttord(Wp,Ws,3,60);
[b,a] = butter(n,Wn,'high');
z = filter(b,a,x);
Z = fft(z,N);
Z2 = (Z/N);
Z1 = Z2(1:N/2+1);%此时选取前半部分，fft之后为对称的双边谱
Z1(2:end-1) = 2*Z1(2:end-1);
figure(4);
subplot(2,1,1);plot(t,z);
title('高通滤波信号（时域）');
subplot(2,1,2);plot(f,Z1);
title('高通滤波信号（频率）');
%带阻滤波
Wp = [0.1,0.5];
Ws = [0.2,0.45]; 
[n,Wn] = buttord(Wp,Ws,3,60);
[b,a] = butter(n,Wn,'stop');
z = filter(b,a,x);
Z = fft(z,N);
Z2 = (Z/N);
Z1 = Z2(1:N/2+1);%此时选取前半部分，fft之后为对称的双边谱
Z1(2:end-1) = 2*Z1(2:end-1);
figure(5);
subplot(2,1,1);plot(t,z);
title('带阻滤波信号（时域）');
subplot(2,1,2);plot(f,Z1);
title('带阻滤波信号（频域）');
%带通滤波
Wp = [0.2,0.3];
Ws = [0.06,0.45];
[n,Wn] = buttord(Wp,Ws,3,60);
[b,a] = butter(n,Wn);
z = filter(b,a,x);
Z = fft(z,N);
Z2 = (Z/N);
Z1 = Z2(1:N/2+1);%此时选取前半部分，fft之后为对称的双边谱
Z1(2:end-1) = 2*Z1(2:end-1);
figure(6);
subplot(2,1,1);plot(t,z);
title('带通滤波信号（时域）');
subplot(2,1,2);plot(f,Z1);
title('带通滤波信号（频域）');

 系统设计

1 MATLAB中的设计

第一步，在matlab中生成x = 5*sin(2*pi*200*t) + 5*sin(2*pi*1000*t) + 5*sin(2*pi*2000*t)的混频信号，然后通过数据处理，得到一组值，方便存放在FPGA的ROM中，实现后面的FPGA滤波。其程序如下：

clear;%清屏
clc;
fs = 8000;N=1024;n=0:N-1;t=n/fs;%sin信号的参数
x = 5*sin(2*pi*200*t) + 5*sin(2*pi*1000*t) + 5*sin(2*pi*2000*t);%混频信号
% c = min(x); %得出信号的最小值，以便知道怎么取正
f = round(x+14)';%根据所得最小值，把信号取正，以便FPGA采用
%c = max(f); %可以得出取正之后数据的最大值，从而知道ROM当中的数据长度
fid = fopen('b.mif','wt');  %将所得的数据存在b.mif当中，（数据文件类型也可以是txt类型）
fprintf(fid,'%g\n',f);
fclose(fid);
[a,b] = size(f); %可以得出数据的位数，好设置ROM的容量
I = a*b;

 

第二步，将低通滤波的参数导出，如下：

%低通滤波
Wp = 200/4000;%0.05
Ws = 400/4000;%0.1
[n,Wn] = buttord(Wp,Ws,0.5,10); %3,60
[b,a] = butter(n,Wn,'low');
  c = min(b);   %将滤波参数取正，以便FPGA当中使用
 
B = round(b*10000+33749)';
A = round(a*10000+33749)';
% c = max(A);
fid = fopen('B.txt','wt');  %将所得的数据存在b.mif当中，（数据文件类型也可以是txt类型）
fprintf(fid,'%g\n',B);
fclose(fid);
 
fid = fopen('A.txt','wt');  %将所得的数据存在b.mif当中，（数据文件类型也可以是txt类型）
fprintf(fid,'%g\n',A);
fclose(fid);

 

2 FPGA(STM32)中的设计

FPGA一共分为４个模块：第一个是pll_clk模块，是用pll IP核进行时钟分频的，用于驱动滤波模块和ROM模块。第二个模块是ROM控制模块，其中ROM的数据来自于matlab中的原始信号的值，对它们进行了取正取整，然后存放在ROM当中，当ROM逐次加一，就可以得到每一位上的值，然后对其直接输出，可以在modelsim看到原始信号的波形，以及输出到滤波模块，进行滤波处理。第三个模块就是滤波模块，通过matlab得到的低通滤波参数进行IIR滤波，最后输出波形。第四个模块，就是顶层模块，将前面三个模块进行元件例化，就可以进行仿真了。（还有一个仿真文件，主要是产生时钟激励）

第一个模块：时钟分频模块

module pll_clk(
    input sys_clk,            //系统时钟
    input sys_rst_n,      //系统复位信号
    output clk_100m,            //输出100mhz时钟
    output rst_n
);
wire locked;
assign rst_n = sys_rst_n & locked;   //当sys_rst_n和locked都为1时，rst_n才为高电平，
                        //所以只有系统稳定时,这个可以用于其他模块的复位状态
my_pll   my_pll_inst(
  .areset(~sys_rst_n),//复位信号和按键复位相连，但锁相环是高电平复位，所以要取反
	.inclk0   (sys_clk)  ,
	.c0       (clk_100m),
	.locked   (locked)         //标志输出稳定
);
endmodule

 第二个模块，rom控制模块，通过地址增加，依次读出数据

//rom的控制模块，用于产生递增的地址信号
module rom_ctrl(
    input clk_100m,
    input rst_n,
    output  q ,        //有效数据输出
    output reg [9:0] addr         //地址输出
);
always @(posedge clk_100m or negedge rst_n)
    begin
        if(!rst_n)
            addr <= 10'd0;   //复位时地址为十进制0
        else
            begin
                if(addr < 1024)
                    addr <= addr + 1; 
                else
                    addr <= 10'd0; 
            end 
    end
//ip核 rom调用
my_rom my_rom_inst(
    .clock(clk_100m),
    .q(q),
    .address(addr)
);   
endmodule

 第三个模块，滤波部分，通过IIR滤波，得到低通信号

module dsp_iir(    //滤波模块
    input clk_100m,
    input rst_n,
    input  [4:0]  q ,        //有效数据输出
    output reg  [24:0]  Q  [4:0]
);
reg [30:0] x_0,x_1,x_2,x_3,x_4,x_5,x_6;
reg [30:0] y_0,y_1,y_2,y_3,y_4,y_5,y_6;
parameter a_0 = 43749;
parameter a_1 = 1;
parameter a_2 = 76873;
parameter a_3 = 9059;
parameter a_4 = 39087;
parameter a = 10000;
parameter b = 33749;
parameter b_0 = 33751;  //1 4  11  14
parameter b_1 = 33755;
parameter b_2 = 33758;
parameter b_3 = 33755;
parameter b_4 = 33751;
always @(posedge clk_100m or negedge rst_n)begin
    if(!rst_n)
        begin
            x_0 <= 0;x_1 <= 0;x_2 <= 0;x_3 <= 0;x_4 <= 0;
            y_0 <= 0;y_1 <= 0;y_2 <= 0;y_3 <= 0;y_4 <= 0;
        end     
    else
        begin
            x_0 = q[0]*a+b;x_1 = q[1]*a+b;x_2 = q[2]*a+b;x_3 = q[3]*a+b;x_4 = q[4]*a+b;//x_5 = q[5];x_6 = q[6];
            y_0 <= (b_0*x_0);
            y_1 <= (b_0*x_1 + b_1*x_0 - a_1*y_0);
            y_2 <= (b_0*x_2 + b_1*x_1 + b_2*x_0 - a_1*y_1 - a_2*y_0);
            y_3 <= (b_0*x_3 + b_1*x_2 + b_2*x_1 + b_3*x_0 - a_1*y_2 - a_2*y_1 - a_3*y_0);            
            y_4 <= (b_0*x_4 + b_1*x_3 + b_2*x_2 + b_3*x_1 + b_4*x_0 - a_1*y_3 - a_2*y_2 - a_3*y_1 - a_4*y_0);
        Q <= y_4;
        end
end 
endmodule

第四个模块，顶层模块：

//系统顶层模块，负责子模块的联级
module dsp_iir_top(
    input sys_clk,
    input sys_rst_n,
    output  Q,
    output  q         //有效数据输出
);
wire [9:0] addr;   //定义地址信号
wire clk_100m;
wire rst_n;
//实例化rom_ctrl
rom_ctrl rom_ctrl(
    .clk_100m(clk_100m),
    .rst_n(rst_n),
    .q(q),
    .addr(addr)
);
pll_clk    u_pll_clk(
    .sys_clk           (sys_clk),
    .sys_rst_n         (sys_rst_n), 
    .rst_n         (rst_n),     
    .clk_100m          (clk_100m       )
    );
dsp_iir    u_dsp_iir(
    .q(q),
    .Q(Q),
    .rst_n         (rst_n),     
    .clk_100m          (clk_100m       )
    );   
endmodule

 第五个模块，测试模块：

`timescale  1ns/1ns
module        dsp_iir_tb      ;
parameter     SYS_PERIOD = 20;  //定义系统时钟周期
reg           sys_clk            ;
reg           sys_rst_n          ;  
wire          clk_100m        ; 
wire          rst_n        ; 
wire          [4:0]  q        ; 
wire          [4:0]  Q        ; 
wire          [9:0]  addr        ; 
always #(SYS_PERIOD/2) sys_clk <= ~sys_clk ;
initial begin
            sys_clk <= 1'b0 ;
            sys_rst_n <= 1'b0 ;
          #(20*SYS_PERIOD)
            sys_rst_n <= 1'b1 ;
          end
//例化ip_pll模块          
rom_ctrl rom_ctrl(
    .clk_100m(clk_100m),
    .rst_n(rst_n),
    .q(q),
    .addr(addr)
);
pll_clk    u_pll_clk(
    .sys_clk           (sys_clk),
    .sys_rst_n         (sys_rst_n), 
    .rst_n         (rst_n),     
    .clk_100m          (clk_100m       )   
    );
dsp_iir    u_dsp_iir(
    .q(q),
    .Q(Q),
    .rst_n         (rst_n),     
    .clk_100m          (clk_100m       )
    );   
//ip核 rom调用
my_rom my_rom_inst(
    .clock(clk_100m),
    .q(q),
    .address(addr)
);
endmodule  

顶层模块图如下：

Modelsim仿真如下：

 

工程文件上传至qq群：868412045