文章目录
前言
提示:FIR(Finite Impulse Response)滤波器,即有限长单位冲激响应滤波器。该滤波器在数字电路的滤波中很常用,本次实验使用Quartus的FIR IP核和MATLAB进行联合仿真,进行FIR IP核的设计与实现。
本篇文章主要介绍如何使用MATLAB和Quartus进行联合仿真,使用MATLAB进行待滤波曲线生成和滤波器参数生成,使用Quartus的FIR IP核设计IP核,并用ROM IP核存储带滤波数据的数据,编写测试文件编写时序并使用Modelsim和MATLAB对结果进行仿真。
本文旨在于进行仿真验证,故不介绍FIR的原理
一、使用MATLAB进行FIR仿真实验
1.1 产生待滤波信号
本次实验需要生成待滤波信号,仿真使用的信号为10M的正弦波,添加50M的正弦波作为噪声, 选择采样率为500M,采样点数1500(该数据可根据需求修改)。
matlab生成上述信号的代码和结果图如下:
%% 初始信号的产生
%正弦信号三要素
A = 60;
f = 10000000; %10M
w = 2*pi*f; %rad/s
p = 0; %rad
A_noise = 0.7*60;
f_noise = 50000000; %50M
w_noise = 2*pi*f_noise; %rad/s
p_noise = 0; %rad
%采样
fs = 500000000; %500MHz %采样频率
d = 1/fs; %s %采样间隔
L = 1500; %信号点数
t = (0:L-1)/fs; %时间轴(X轴)数据
s1=A*sin(w*t+p); %正弦信号
s2=A_noise*sin(w_noise*t+p_noise); %正弦信号
s = s1 + s2; %融合后的信号
%用均值为0,方差为4的白噪声扰乱该信号
%x = s + 2*randn(size(t));
figure(1);
plot(t(1:1000),s(1:1000)) %因为频率较高,只画一部分
%可以看出很难确定频率分量
xlabel('时间/s');
ylabel('幅度');
使用MATLAB自带的FFT函数对其进行验证,程序和结果如下。从结果图中可以看出原始信号经过FFT后显示有一个10M的信号和一个50M的信号
%% FFT
%计算双侧频谱 P2。然后基于 P2 和偶数信号长度 L 计算单侧频谱 P1。
y = fft(s);
P2 = abs(y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);
%定义频域 f 并绘制单侧幅值频谱 P1。
f = fs*(0:L/2)/L;
%f = Fs*(0:(L/2))/L;
figure(2);
plot(f,P1);
title('Magnitude');
xlabel('f (Hz)')
ylabel('|P1(f)|')
1.2 使用MATLAB自带滤波器工具箱设计滤波器
在MATLAB下方的 命令行窗口 中直接输入filterDesigner并回车,即可打开MATLAB自带的滤波器设计工具箱,该工具箱中的操作步骤如下:
步骤一:选择其为定点8位,具体位数取决于你需要处理的数据位数,本次设计数据在-128~127之间,故此处选择的位数为8位,具体设置如下图所示。
这部分需要与后续的ROM IP核以及FIR IP核中的一些参数相同
步骤二:选择低通滤波器,阶数Specify order为32位(阶数越高滤波效果越好),频率选择处Fs为采样率,Fpass为低通滤波时开始转折下降的频率,Fstop为第一次下降到底的频率。具体操作如下图所示
步骤三:导出FIR滤波器系数。后续FIR仿真和FIR IP核中均需要滤波器系数,故需要将其导出。导出时会默认将滤波器系数导出为变量Num,如果在Workspace中已经有此名称的数据,则需要更换名称,否则会报错。导出方法如下如所示。
1.3 MATLAB仿真结果
在滤波器系数导出后,直接将原曲线与Num(滤波器系数)直接进行卷积即可。此部分程序如下
%% FIR
%先使用filterDesigner导出滤波器系数Num
out = conv(s,Num); %卷积
figure(4)
subplot(2,1,1)
plot(t(1:1000),s(1:1000)) %因为频率较高,只画一部分
title('滤波前');
xlabel('时间/s');
ylabel('幅度');
subplot(2,1,2);
plot(t(1:1000),out(1:1000)) %因为频率较高,只画一部分
title('滤波后');
xlabel('时间/s');
ylabel('幅度');
利用和1.1节相同方法对滤波后的曲线进行FFT,得到的曲线如下图(此部分代码和之前非常类似,故不再另外放置代码)
因为卷积会导致整个数据增多,故fft时的数据只使用(1:1500)即可,即为y = fft(out(1:L));
二 基于FIR IP核的设计与实现
2.1 使用ROM IP核生成待滤波信号
鉴于待测信号难以使用代码生成,故选择将数据导入到ROM IP核中生成。ROM IP核支持hex文件和mif文件,本次使用的是mif文件。直接使用MATLAB生成mif文件,这部分参考小梅哥的生成mif文件程序。
因为mif文件中的数据为hex,不支持浮点数和负数。上一章中使用的数据为double型且有正数和负数,故需要先对数据进行修改再写入为mif文件。
先将数据的最小值找出,再让原曲线每个数据减去该值(因为该值为负数),最后使用round函数进行取整,具体程序如下
因为之前设置的位数为8,所以要保证这部分的数据在0~255范围内
%% 将s数据全变为正数
min_s = min(s);
s_mif = zeros(size(s)) ;
for i = 1:L
s_mif(i) = round(s(i) - min_s); %对小数四舍五入以取整
if s_mif(i) <0 %强制将负1置0,
s_mif(i) = 0;
end
end
figure(3)
plot(t(1:1000),s_mif(1:1000)) %因为频率较高,只画一部分
%可以看出很难确定频率分量
xlabel('时间/s');
ylabel('幅度');
%% 产生Quartus生成ROM的IP核所需的.mif文件
%生成数据为s
fild = fopen('D:\FPGA_STUDY\program\my_pro\sin10M_50M_500MSPS.mif','wt');%创建mif文件
%写入mif文件文件头
fprintf(fild, '%s\n','WIDTH=8;');%位宽
fprintf(fild, '%s\n\n','DEPTH=1500;');%深度
fprintf(fild, '%s\n','ADDRESS_RADIX=UNS;');%地址格式
fprintf(fild, '%s\n\n','DATA_RADIX=HEX;');%数据格式
fprintf(fild, '%s\t','CONTENT');%地址
fprintf(fild, '%s\n','BEGIN');%
str = round(s);
for i = 1:L
% addr : data;
fprintf(fild, '\t%g\t',i-1);%地址,从0开始编码
fprintf(fild, '%s\t',':');
fprintf(fild, '%x',s_mif(i));
%fprintf(fild, '%d',str(i));
fprintf(fild, '%s\n',';');
end
fprintf(fild, '%s\n','END;');
2.2 ROM IP核使用说明
本次实验使用的是单端口ROM,ROM位宽为8,深度选2048。导入2.1节生成的mif文件。ROM IP核的单端口ROM操作非常简单,只需要不断的增加地址即可。实验在程序中添加了rden信号,但是一直给高,并没有起作用。
2.3 FIR IP核使用说明
本次实验使用的FIR IP核为FIR Compiler v13.1,具体设置步骤如下。
注意:FIR IP核最好在工程目录下生成,否则进行仿真时可能会报错。
点击step1,进入下面界面。点击Edit Coefficient Set进入滤波器参数设置界面
此步骤即将从MATLAB的滤波器设计工具箱中得到的滤波器参数Num导入到IP核中。将MATLAB的Num导出为TXT文件,每个参数占据一行,此处参数无需修改为正数整数之类,直接用导出的参数即可。
此部分设置数据位宽,这部分需要和之前在MATLAB中的滤波器工具箱中设置的参数一致。数据类型选择无符号数,因为本次我们导出的.mif文件中的数据是无符号数据。
可以看到下面的提示中说数据的输入和输出均为一个时钟周期,这部分在后续导出数据中有用。
step2的内容就直接全部勾选,之后生成。生成的时候可能需要不断的生成,即一次生成可能会卡在某个位置不动,此时直接重复step3即可。
FIR IP核生成内容可参考正点原子的教程,此处链接:正点原子FIR IP核教程
2.4 仿真测试文件编写
仿真测试文件由两个部分组成,即顶层文件核仿真文件,顶层文件主要设计时序,仿真测试文件主要用于导出数据。顶层文件为fir_top.v,仿真测试文件(testbench)为fir_tb.v
FIR IP核的端口和时序如下图
ast_sink_error需要给2‘b0,否则在仿真时ast_source_valid信号无法拉高
FIR IP核单通道时序如下图,实验中将ast_sink_valid信号一直拉高,ast_sink_ready置一时开始控制ROM IP核地址线自加,即开始向FIR IP核中写入待滤波数据。
module fir_top(
input Clk, //50M
input Rst_n,
output [17:0] fir_data, //经过FIR之后的数据
output [7:0] data_sin,
output ast_source_valid
);
reg [10:0] addr;
//ROM例化
rom rom_u(
.address(addr),
.clock(Clk),
.q(data_sin),
.rden(1'b1) //当FIR模块为高时,使能rden,此后开始输出数据
//rden信号拉低后,下一个时候上升沿检测到,但这次会取出一个数据,在拉低之后的第二个时钟上升沿会停止取数,数据线上的值会保持不变,而不是清零
);
wire ast_sink_valid;
wire ast_sink_ready; //当ast_sink_ready信号为高时表示FIR模块准备完成,此时可向FIR模块中写入数据
//wire ast_source_valid;
reg flag_ready;
assign ast_sink_valid = (addr == 0) ? 1'b0 : 1'b1; //当rom的地址为0时ast_sink_valid信号为0,否则为1
wire [1:0] source_error;
//fir例化
fir fir_u(
.clk(Clk),
.reset_n(Rst_n),
.ast_sink_data(data_sin),
.ast_sink_valid(ast_sink_valid), //input
// Asserted when input data is valid. When ast_sink_valid is not asserted, the FIR processing is stopped if new data is required and no data is left in the AvalonST input FIFO. Otherwise, the FIR processing continues.
.ast_source_ready(1'b1), //input
//Asserted by the downstream module if it is able to accept data
.ast_sink_error(2'd0), //一定要加,否则ast_source_valid信号无法拉高
.ast_source_data(fir_data),
.ast_sink_ready(ast_sink_ready), //output
//Asserted by the FIR filter when it is able to accept data in the current clock cycle
.ast_source_valid(ast_source_valid), //output
//Asserted by the FIR filter when there is valid data to output
.ast_source_error()
);
always@(posedge Clk or negedge Rst_n)begin
if(!Rst_n)begin
addr <= 11'd0;
//flag_ready <= 1'b0;
end
else if(ast_sink_ready == 1)
if(addr >= 1499)begin
addr <= 1;
//flag_ready <= 1'b0;
end
else begin
//flag_ready <= 1'b1;
addr <= addr + 1'b1;
end
else
addr <= addr;
end
//
// reg [15:0] cnt;
// always@(posedge Clk or negedge Rst_n)begin
// if(!rst_n)
// end
endmodule
仿真时需要将数据导出为txt,以便将数据导出到MATLAB中进行FFT分析,程序如下。在ast_source_valid信号拉高后开始写入数据,因为前文说过每一个周期读出一个有效数据,故需要每次延时一个周期后向文件中写入数据。
`timescale 1ns/1ns
`define clk_period 20
module fir_tb;
reg clk;
reg Rst_n;
wire [17:0] fir_data;
wire [7:0] data_sin;
fir_top fir_top_u(
.Clk(clk), //50M
.Rst_n(Rst_n),
.fir_data(fir_data), //经过FIR之后的数据
.data_sin(data_sin),
.ast_source_valid(ast_source_valid)
);
initial clk = 1;
always #(`clk_period/2) clk = ~clk;
integer i;
integer w_file;
initial begin
Rst_n = 0;
#101;
Rst_n = 1;
wait(ast_source_valid == 1); //等待ast_source_valid信号拉高
w_file = $fopen("fir_data.txt");
for(i=0;i<1500;i=i+1)
begin
wait(clk);
#(`clk_period);
$fdisplay(w_file,"%d",fir_data);
if(i == 1500)
$fclose(w_file);
end
#(`clk_period * 2000);
$stop;
end
endmodule
2.5 仿真结果
在进行仿真时需要在设置中导入fir_tb.v文件和待滤波数据.mif文件,其设置如如下。
先使用Modelsim观察波形,选择fir_data和data_sin,将两者均设置为unsigned,并设置为模拟显示,如下图。
点击之后的波形如图所示,可以看到滤波器输出的波形确实已经为正弦波了。
将导出的fir_data导出到MATLAB中并绘图和FFT,可以看到和直接使用MATLAB的卷积运算滤波效果差不多,FFT结果也一致,但是幅值大了很多,这部分怀疑是由于将正弦波加了一个直流分量导致的,这部分后续可以使用一定占空比方波进行验证。
将TXT中的数据导入MATLAB如下图所示
该部分数据为直接导出的数据,并未去除直流分量
对去除直流分量的数据进行FFT,显示其频率为10MHz,即仿真验证成功。
最后一部分MATLAB程序如下,此部分去掉了数据的直流分量。
%% 测试FIR IP核仿真结果
fir_data = fir_data2(50:1499);
L_fir = length(fir_data);
fir_data_mean = fir_data - mean(fir_data); %去除数据直流分量
%计算双侧频谱 P2。然后基于 P2 和偶数信号长度 L 计算单侧频谱 P1。
y = fft(fir_data_mean);
P2 = abs(y/L_fir);
P1 = P2(1:L_fir/2+1);
P1(2:end-1) = 2*P1(2:end-1);
%定义频域 f 并绘制单侧幅值频谱 P1。
f = fs*(0:L_fir/2)/L_fir;
figure(6);
plot(t(1:L_fir),fir_data) %因为频率较高,只画一部分
title('滤波后');
xlabel('时间/s');
ylabel('幅度');
%f = Fs*(0:(L/2))/L;
figure(7);
plot(f,P1);
title('Magnitude');
xlabel('f (Hz)')
ylabel('|P1(f)|')
20211224添加
之前说过FIR 的IP核直接输出的波形的幅值会增加,最近研究发现有FIR IP核设计中有一个选项可以改变输出位数,进而将输出信号的幅值修改到目标位数,以满足后续实验需要。
第一步将Output Number System调节为Custom Resolution
第二步将保留位数修改为目标位数,此处我选择的是8位。当然如果按上面看是保留了高8位,我猜测如果直接输出是输出高8位输出也是可以的,这部分还有待验证
总结
整个仿真流程如上图,此次仿真流程只是验证了FIR IP核可以进行低通滤波。还有很多直接编写Verilog程序实现低通滤波的,并没有使用FIR IP核,这部分后续也会研究研究。如果上述过程存在疑问欢迎留言讨论,转载请注明出处。