FPGA设计实现FIR滤波器仿真验证
要求：设计一个15阶(长度为16)的低通线性相位FIR滤波器，采用blackman窗，截止频率500Hz,采样频率2KHz,采用FPGA全串行结构滤波器实现，滤波器系数量化位宽12bit，输入数据位宽12bit，输出数据位宽29bit,系统时钟16KHz。(数据速率刚好2KHz)
目录
1.首先采用Matlab求出符合要求的滤波器的系数，并观察12bit量化后的滤波器幅频响应曲线（代码如下）。
2.采用Matlab仿真产生待测试数据
3.编写Verilog代码实现FPGA设计滤波器结构
4.编写仿真测试文件
5.将仿真滤波器滤波后输出的数据写入txt文件中。
6.采用Matlab分析经过滤波后的数据，分析滤波性能。
7.总结
1.首先采用Matlab求出符合要求的滤波器的系数，并观察12bit量化后的滤波器幅频响应曲线（代码如下）。
%fir8serial.m
function hn = E4_7_Fir8Serial
N = 16;%滤波器长度
fs = 2000;%滤波器采样频率
fc = 500;%低通滤波器的截止频率
B = 12;%量化位数%生成窗函数
w_black = blackman(N)';
%采用fir1函数设计FIR滤波器
b_black = fir1(N-1,fc*2/fs,w_black);
%量化滤波器系数
Q_black = round(b_black/max(abs(b_black))*(2^(B-1)-1));
hn = Q_black;
%转换成16进制数补码
Q_h = dec2hex(Q_black+2^B*(Q_black<0));%求滤波器的幅频响应
m_black = 20*log10(abs(fft(b_black,1024)));m_black = m_black - max(m_black);%将最大值归位0
Ql_black = 20*log10(abs(fft(Q_black,1024)));Ql_black = Ql_black - max(Ql_black);%设置幅频响应横坐标单位为Hz
x_f = [0:(fs/length(m_black)):fs/2];
m1 = m_black(1:length(x_f));
m2 = Ql_black(1:length(x_f));%绘制幅频响应曲线
plot(x_f,m1,'-',x_f,m2,'-.');
xlabel('频率(Hz)');ylabel('幅度(dB)');
legend('未量化','12bit量化');
grid;
运行结果为
可以发现12bit量化符合要求
滤波器系数如下
2.采用Matlab仿真产生待测试数据
并通过设计的滤波器验证理论上的滤波效果同时将生成的测试数据写入txt文件中供后面Modelsim仿真时读取
代码如下：
%E4_7_NoiseAndCarrier.m
f1 = 200 ;%信号1频率200Hz
f2 = 800 ;%信号2频率800Hz
Fs = 2000;%采样频率2KHz
N = 12 ;%量化数据位宽12bit 输入数据12bit%产生信号
t = 0:1/Fs:1;
c1 = 2*pi*f1*t;
c2 = 2*pi*f2*t;
s1 = sin(c1);%产生正弦信号
s2 = sin(c2);%产生正弦信号
s = s1 + s2;%产生两个单载波合成后的信号
%产生随机序列信号
noise=randn(1,length(t));%产生高斯白噪声序列%归一化处理
s = s / max(abs(s));
noise = noise / max(abs(noise));%12bit量化
Q_noise = round(noise*(2^(N-1)-1));
Q_s = round(s*(2^(N-1)-1));%调用自己设计的滤波器函数对信号进行滤波
hn = fir8serial;
Filter_noise = filter(hn,1,Q_noise);
Filter_s = filter(hn,1,Q_s);%求信号的幅频响应
%滤波之前
m_noise = 20*log10(abs(fft(Q_noise,1024))); m_noise = m_noise - max(m_noise);
m_s = 20*log10(abs(fft(Q_s,1024))); m_s = m_s - max(m_s);
%滤波之后
Fm_noise = 20*log10(abs(fft(Filter_noise,1024)));Fm_noise = Fm_noise - max(Fm_noise);
Fm_s = 20*log10(abs(fft(Filter_s,1024)));Fm_s = Fm_s - max(Fm_s);
%滤波器本身的幅频响应
m_hn = 20*log10(abs(fft(hn,1024)));m_hn = m_hn - max(m_hn);%设置幅频响应的横坐标为Hz
x_f = [0:(Fs/length(m_s)):Fs/2];
%只显示正频率部分
mf_noise = m_noise(1:length(x_f));
mf_s = m_s(1:length(x_f));
Fmf_noise = Fm_noise(1:length(x_f));
Fmf_s = Fm_s(1:length(x_f));
Fm_hn = m_hn(1:length(x_f));%绘制幅频响应曲线
subplot(211)
plot(x_f,mf_noise,'-.',x_f,Fmf_noise,'-',x_f,Fm_hn,'--');
xlabel('频率(Hz)');ylabel('幅度(dB)');title('Matlab仿真白噪声滤波前后的频谱');
legend('输入信号频谱','输出信号频谱','滤波器响应');
grid;subplot(212);
plot(x_f,mf_s,'-.',x_f,Fmf_s,'-',x_f,Fm_hn,'--');
xlabel('频率(Hz)');ylabel('幅度(dB)');title('Matlab仿真合成单频信号滤波前后频谱');
legend('输入信号','输出信号','滤波器响应');
grid;%将生成的数据以10进制格式写入txt文档
fid = fopen('D:\test\fir\4_7\E4_7_Int_noise.txt','w');
fprintf(fid,'%8d\r\n',Q_noise);
fprintf(fid,';');
fclose(fid);fid = fopen('D:\test\fir\4_7\E4_7_Int_s.txt','w');
fprintf(fid,'%8d\r\n',Q_s);
fprintf(fid,';');
fclose(fid);%将生成的数据以二进制形式写入txt文档
fid = fopen('D:\test\fir\4_7\E4_7_Bin_noise.txt','w');
for i = 1:length(Q_noise)B_noise = dec2bin(Q_noise(i)+(Q_noise(i)<0)*2^N,N);%N为位数for j=1:Nif B_noise(j)=='1'tb = 1;elsetb = 0;end fprintf(fid,'%d',tb);end fprintf(fid,'\r\n');
end
fprintf(fid,';');
fclose(fid);fid = fopen('D:\test\fir\4_7\E4_7_Bin_s.txt','w');
for i=1:length(Q_s)B_s = dec2bin(Q_s(i)+(Q_s(i)<0)*2^N,N);%如果为负数则转换为它所对应的正数for j=1:Nif B_s(j)=='1'tb = 1;else tb = 0;end fprintf(fid,'%d',tb);end fprintf(fid,'\r\n');
end
fprintf(fid,';');
fclose(fid);运行结果：
我们可以发现当输入信号为白噪声时输出800Hz以后的频率衰减已经达到-40dB了
当输入信号为200Hz和800Hz的叠加信号时，滤波后的输出只剩200Hz的信号了，800Hz被滤除了。
3.编写Verilog代码实现FPGA设计滤波器结构
/*
*
*@Author: X-Z
*@Date:2023-09-17 00:45:22
*@Function:全串行FIR数字滤波器 加法器位宽13bit,量化12bit扩展一位 乘法器输入位宽12bit和13bit输出数据位宽25bit
*//*设计一个15阶(长度2为16)的低通线性相位的FIR滤波器，采用布莱克曼窗，截止频率500Hz,采样频率2KHZ,采用FPGA实现全串行结构的滤波器，系数量化位数12Bit,输入数据位宽12bit，输出数据位宽29bit系统时钟16KHz
*/module FirFullSerial(input wire clk ,//FPGA系统时钟频率16KHzinput wire rst_n ,input wire signed [11:0] Xin ,//数据输入频率2KHzoutput wire signed [28:0] Yout //滤波器输出数据
);//例化有符号数乘法器IP核reg signed [11:0] coe ;//滤波器为12bit量化数据wire signed [12:0] add_s;//输入为12bit量化滤波器系数，两个对称项系数相加需要13bit位宽wire signed [24:0] Mout ;
mul u_mult(.aclr (~rst_n),.clock (clk ),.dataa (coe ),.datab (add_s ),.result (Mout ));//例化有符号数加法器IP核,对输入数据进行1位数据位的扩展，输出结果为13bit数据reg signed [12:0] add_a;reg signed [12:0] add_b;
add u_adder(.dataa (add_a),.datab (add_b),.result (add_s));//3位计数器，计数周期为8，刚好为数据输入的速率2KHZreg [2:0] cnt;always @(posedge clk or negedge rst_n)beginif(!rst_n)cnt <= 3'd0;elsecnt <= cnt + 1'd1; end//将数据移入移位寄存器Xin_Reg中reg [11:0] Xin_Reg [15:0];//16个12bit位宽的存储器，用来存储滤波器系数reg [3:0] i,j;always @(posedge clk or negedge rst_n)beginif(!rst_n)begin//初始化寄存器为0for(i=0;i<15;i=i+1)Xin_Reg[i] <= 12'd0;endelse beginif(cnt==3'd7)beginfor(j=0;j<15;j=j+1)Xin_Reg[j+1] <= Xin_Reg[j];Xin_Reg[0] <= Xin; endendend//将对称系数的输入数据相加，同时将对应滤波器的系数送人乘法器
//需要注意的是以下程序只是用了一个加法器和一个乘法器资源
//以8倍速率调用乘法器IP核，由于滤波器长度为16，系数具有对称性，所以
//可以在一个时钟周期内完成所有8个滤波器系数与对应数据的乘法运算
//为保证加法运算不溢出，输入输出数据均扩展为13bit
//乘法器使用了1级流水线，加法器没有使用流水线会立即输出相加结果
/*12bit量化滤波器系数0 -3 15 46 -117 -263 590 2047 2047 590 -263 -117 46 15 -3 0
*/always @(posedge clk or negedge rst_n)beginif(!rst_n)beginadd_a <= 13'd0;add_b <= 13'd0;coe <= 12'd0;//滤波器系数endelse begincase(cnt)3'd0 : beginadd_a <= {Xin_Reg[0][11],Xin_Reg[0]};//将数据扩展一位符号位add_b <= {Xin_Reg[15][11],Xin_Reg[15]};coe <= 12'd0;//c0end3'd1 : beginadd_a <= {Xin_Reg[1][11],Xin_Reg[1]};add_b <= {Xin_Reg[14][11],Xin_Reg[14]};coe <= -12'd3;//c1end3'd2 : beginadd_a <= {Xin_Reg[2][11],Xin_Reg[2]};add_b <= {Xin_Reg[13][11],Xin_Reg[13]};coe <= 12'd15;//c2end3'd3 : beginadd_a <= {Xin_Reg[3][11],Xin_Reg[3]};add_b <= {Xin_Reg[12][11],Xin_Reg[12]};coe <= 12'd46;//c3end3'd4 : beginadd_a <= {Xin_Reg[4][11],Xin_Reg[4]};add_b <= {Xin_Reg[11][11],Xin_Reg[11]};coe <= -12'd177;//c4end3'd5 : beginadd_a <= {Xin_Reg[5][11],Xin_Reg[5]};add_b <= {Xin_Reg[10][11],Xin_Reg[10]};coe <= -12'd263;//c5end3'd6 : beginadd_a <= {Xin_Reg[6][11],Xin_Reg[6]};add_b <= {Xin_Reg[9][11],Xin_Reg[9]};coe <= 12'd590;//c6end 3'd7 : beginadd_a <= {Xin_Reg[7][11],Xin_Reg[7]};add_b <= {Xin_Reg[8][11],Xin_Reg[8]};coe <= 12'd2047;//c7enddefault : beginadd_a <= 13'd0;add_b <= 13'd0;coe <= 12'd0;//滤波器系数endendcaseendend//对滤波器系数与输入数据的乘法结果进行累加，并输出滤波后的数据
//考虑到乘法器及累加器的延时，需要计数器为2时对累加器清零，同时输出
//滤波器结果数据，类似的时延长度一方面可以通过精确计算获得，但更好的方法是
//通过行为仿真查看reg signed [28:0] sum ;reg signed [28:0] yout; always @(posedge clk or negedge rst_n)beginif(!rst_n)beginsum <= 29'd0;yout <= 29'd0;endelse beginif(cnt==3'd2)beginyout <= sum ;//同时cnt==1时最后一个数运算完成但cnt==2时才有效；输出yout有效时cnt==3sum <= Mout;//ct==2时输出第一个乘法器运算结果endelse sum <= sum + Mout; //经过1级流水线，cnt==1时全部加完cnt==2时输出累加结果有效endendassign Yout = yout;endmodule4.编写仿真测试文件
`timescale 1ns/1ps
module FirFullSerial_tb();reg clk ;reg rst_n ;reg [11:0] Xin ;wire clk_data;//数据时钟，速率为系统时钟clk速率的1/8wire [28:0] Yout ;//例化滤波器模块
FirFullSerial u_FirFullSerial(/*input wire */ .clk (clk ) ,//FPGA系统时钟频率16KHz/*input wire */ .rst_n(rst_n) ,/*input wire signed [11:0] */ .Xin (Xin ) ,//数据输入频率2KHz/*output wire signed [28:0] */ .Yout (Yout ) //滤波器输出数据
);parameter clk_period = 62500;//设置时钟信号周期(频率)：16KHzparameter clk_period_data = clk_period * 8;parameter clk_half_period = clk_period / 2;parameter clk_half_period_data = clk_period_data / 2;parameter data_num = 2000;//仿真数据长度parameter time_sim = data_num * clk_period_data;//仿真时间initial beginclk <= 1'b1;rst_n = 1'b0;#(5*clk_period)rst_n = 1'b1;//设置仿真时间#time_sim $finish;Xin = 12'd10;//输入信号赋初值end//产生时钟信号always #clk_half_period clk = ~clk;//16KHzreg [2:0] cn_clk = 0;always @(posedge clk) cn_clk <= cn_clk + 3'd1;assign clk_data = cn_clk[2];//时钟周期与系统时钟对齐；2KHz//从外部txt文件(SinIn.txt)中读入数据作为测试数据integer pattern;reg [11:0] stimuls [1:data_num];initial begin//文件必须放在"工程目录\simulation\modelsim"路径下$readmemb("E4_7_Bin_noise.txt",stimuls);//$readmemb("E4_7_Bin_s.txt",stimuls);pattern = 0;repeat(data_num)begin@(posedge clk_data);pattern = pattern + 1;Xin = stimuls[pattern];endend//将仿真数据dout写入外部txt文件中(out.txt)integer file_out;initial begin//放在"工程目录\simulation\modelsim"路径下file_out = $fopen("E4_7_Noiseout.txt");//file_out = $fopen("E4_7_sout.txt");if(!file_out)begin$display("could not open file!");$finish;endendwire rst_write;wire signed [28:0] dout_s;assign dout_s = Yout;//将dout转换成有符号数据assign rst_write = clk_data & rst_n;//产生写入时钟信号，复位状态下不写入数据always @(posedge rst_write)$fdisplay(file_out,"%d",dout_s);endmodule 5.将仿真滤波器滤波后输出的数据写入txt文件中。
Modelsim仿真结果：
通过仿真我们可以发现大致实现了滤波性能，但还是不够直观，因此我们采用Matlab进行分析滤波性能。
6.采用Matlab分析经过滤波后的数据，分析滤波性能。
代码如下：
%E4_7_NoiseAndCarrierOut.m
f1 = 200 ;%信号1频率200Hz
f2 = 800 ;%信号2频率800Hz
Fs = 2000;%采样频率2KHz
N = 12 ;%量化位数%从文本文件中读取数据
%测试输入数据分别放在Noise_in和S_in变量中
fid = fopen('D:\test\fir\4_7\E4_7_Int_Noise.txt','r');
[Noise_in,N_n] = fscanf(fid,'%lg',inf);
fclose(fid);fid = fopen('D:\test\fir\4_7\E4_7_Int_s.txt','r');
[S_in,N_n] = fscanf(fid,'%lg',inf);
fclose(fid);%滤波后的输出结果数据放在Noise_out和S_out变量中
fid = fopen('D:\test\fir\4_7\FirFullSerial\prj\simulation\modelsim\E4_7_Noiseout.txt','r');
[Noise_out,N_count] = fscanf(fid,'%lg',inf);
fclose(fid);fid = fopen('D:\test\fir\4_7\FirFullSerial\prj\simulation\modelsim\E4_7_sout.txt','r');
[S_out,N_count] = fscanf(fid,'%lg',inf);%以浮点数形式全部读取文件中的数据
fclose(fid);%归一化处理
Noise_out = Noise_out / max(abs(Noise_out));
S_out = S_out / max(abs(S_out));
Noise_in = Noise_in / max(abs(Noise_in));
S_in = S_in /max(abs(S_in));%求信号的幅频响应
out_noise = 20*log10(abs(fft(Noise_out,1024)));out_noise = out_noise - max(out_noise);
out_s = 20*log10(abs(fft(S_out,1024)));out_s = out_s - max(out_s);
in_noise = 20*log10(abs(fft(Noise_in,1024)));in_noise = in_noise - max(in_noise);
in_s = 20*log10(abs(fft(S_in,1024)));in_s = in_s - max(in_s);%滤波器本身的幅频响应
hn = fir8serial;
m_hn = 20*log10(abs(fft(hn,1024)));m_hn = m_hn - max(m_hn);%设置幅频响应横坐标单位为Hz
x_f = [1:(Fs/length(out_noise)):Fs/2];
%只显示正频率部分的幅频响应
mf_out_noise = out_noise(1:length(x_f));
mf_out_s = out_s(1:length(x_f));
mf_in_noise = in_noise(1:length(x_f));
mf_in_s = in_s(1:length(x_f));
mf_hn = m_hn(1:length(x_f));%绘制幅频响应特性曲线
figure(1);
subplot(211);
plot(x_f,mf_in_noise,'--',x_f,mf_out_noise,'-',x_f,mf_hn,'-.');
xlabel('频率(Hz)');ylabel('幅度(dB)');
title('FPGA仿真白噪声滤波前后幅频响应曲线');
legend('输入信号频谱','输出信号频谱','滤波器响应');
grid;
subplot(212);
plot(x_f,mf_in_s,'--',x_f,mf_out_s,'-',x_f,mf_hn,'-.');
xlabel('频率(Hz)');ylabel('幅度(dB)');title('FPGA仿真单频叠加信号滤波前后幅频响应曲线');
axis([0 1000 -100 0]);
legend('输入信号幅频响应','输出信号频谱','滤波器响应');
grid;%绘制时域波形
%设置显示数据的范围
t = 0:1/Fs:50/Fs;t=t*1000;
%length(t)
t_in_noise = Noise_in(1:length(t));
t_in_s = S_in(1:length(t));
t_out_s = S_out(1:length(t));
t_out_noise = Noise_out(1:length(t));figure(2);
subplot(211);
%length(hn);
plot(t,t_in_noise,'--',t,t_out_noise,'-');
xlabel('时间(ms)');ylabel('幅值');title('FPGA仿真白噪声滤波前后的时域波形');
legend('输入信号波形','输出信号波形');
grid;
subplot(212);
plot(t,t_in_s,'--',t,t_out_s,'-');
xlabel('时间(ms)');ylabel('幅值');
title('FPGA仿真单频叠加信号滤波前后时域波形');
legend('输入信号波形','输出信号波形');
grid;运行结果为：
我们发现采用FPGA滤波前后的波形和Matlab仿真波形图基本一致，唯一不同的是Matlab的滤波效果更好一些，FPGA衰减没有Matlab仿真的大，这是因为FPGA中存在字长效应以及数据处理时存在截尾效应。因此可以发现理论与实际的差距，但最终设计满足要求。
注：自学完杜勇的FIR滤波器所写，能力有限，不喜勿喷。
本文来自互联网用户投稿，文章观点仅代表作者本人，不代表本站立场，不承担相关法律责任。如若转载，请注明出处。 如若内容造成侵权/违法违规/事实不符，请点击【内容举报】进行投诉反馈！