MATLAB设计FPGA实现联合ISE和Modelsim仿真的FIR滤波器设计

1、MATLABS计，FPGA实现，联合ISE和Modelsim仿真的FIR滤波器设计摘要：本文首先利用 MATLAE产生两个频率不一样的正弦信号，并将这两个正弦信号相加，得到一个混叠的波形；然后利用 MATLAB计一个FIR低通滤波器，并由 Verilog 实现，联合ISE和Modelsim仿真，实现滤除频率较高的信号，并将滤波后的数据送到MATLAB分析。本文首先介绍FIR滤波器的基本原理，然后从目前最流行的的设计方式入手，逐步介绍FIR滤波器的设计，其流程图如图1所示：图1FIR滤波器的基本原理：数字滤波器是从分析信号中提取用户需要的信息，滤去不需要的信号成分或者干扰成分。数字滤波器的数学模

2、型：线性时不变数字滤波器的数学模型在时域中可以用线性常 系数差分方程给出：NMyndkLyn-k Pq_xn-qk=1q 二0M其等效的Z变换为:PqLzT1dkLz*k d当dk全为0时，Z域系统函数只有零点，数字滤波器的单位脉冲响应有限，通常这用滤波器被称为有限冲激响应( FIR)滤波器。它的输出只与当前的输入和过去的输 入有关，跟过去的输出无关。FIR滤波器的一个重要特性是具有线性相位，即系统的相 异和频率成正比，可达到无失真传输。F1R滤波器可以用方框图的形式方便地表示.用方框團法衷示有戏下几个好处：可通过 观察法容易地写出算法.通过调黑框區得到不同算法的饕效框图可以容易地确定硬件的需

3、 求；此外.还可以从传输函敷所土成的框图宜援御到名种等效表寿*其亘接形武站方框图如图63所示。”曲；1汀f卩制如(01 M2)u曲卅_即仃叫 4 X耐图&3亘接形式的FIR滤波器利用转置定理，可汉将直接型转化为其等效的转置式鬥R滋波嚣.挂结构如图&亠表不* 转SHR的优点在于不需要给点町提供额外的移位寄存器，而且也没有必要为达到高速牡理 给乘积的加法誥脈加额外的流水线f图6斗 转置式的FIR滤波鶉更加详细的介绍FIR滤波器的理论算法知识，请参考别的先关书籍资料、二、 FIR滤波器的设计及验证：1)由MATLAB生FIR滤波器的输入数据，并存放到.txt文件中用MATLAB生两个频率分别为：0.

4、5MHz和2MHz的正弦信号，并将这两个正弦信号叠加，得到一个混叠的信号，并将这个混叠信号的波形数据存放到C盘根目录下面的matlab_wave_data.txt 文件中。MATLAB 现代码如下：clear all;fs = 25000000; % 25M 采样率t = 0:1/fs:0.0002; % 共 0.0002 * 25000000 = 5000 个点fl = 500000;f2 = 2000000;sig nail = si n(2*pi*f1*t); %频率为 0.5MHz 的正弦信号%由于正、余弦波形的值在0,1之间，需要量化成16bit,先将数值放大y2 = fix(163

5、84 + (2A14 - 1) * sig nail); %量化到 15bity3 = fix(16384 + (2人14 - 1) * sig nal2); %量化到 15bity1 = y2 + y3; % 量化到 16bit%再将放大的浮点值量化，并写道到存放在C盘的文本中fid = fope n(c:/matlab_wave_data.txt, wt);fprintf(fid, %xn, y1); %以16进制的格式将数据存放到.txt 文件中。如果不用 16进制，后面%读数会出错。fclose(fid);plot(y2, b);hold on;plot(y3, b);hold on;

6、plot(y1, r);legend(0.5MHz 正弦,2MHz 正弦,两者叠加);title(MATLAB产生的两个正弦信号的叠加波形);存放数据的文件路径、名称已经文本的内容如图3所示:图3用MATLAB读C盘根目录下的 matlab_wave_data.txt 文件，验证存入的波形数据是否 正确，MATLAB代码如下：fid = fope n(c:/matlab_wave_data.txt,r);for i = 1 : 5001;% 共有 5001 个数据num(i) = fscanf(fid, %x, 1);%从fid所指的文件中，以16进制的方式读出一个数据endfclose(fi

7、d);figure(2);plot( nu m,b);legend(MATLAB从txt文件中读出的原始叠加波形数据);title( 直接回读MATLAB生的两个正弦信号的叠加波形);显示的波形如图4所示:图4对比图4和图2中的叠加波形，可知以上操作的正确性。2) 用MATLAB计FIR滤波器输入信号是频率别为 0.5MHz和2MHz的正弦信号的叠加，我们的任务是设计一个低通滤 波器滤除掉2MHz的干扰信号。因此，我们可以设计一个采样率为25MHz的低通滤波器，其通带带宽为1MHz阻带宽度为2MHz通带内纹波抖动为 1dB，阻带下降为 80dB。在MATLAB勺命令窗口输入：fdatool命令

8、并回车，打开 FDATool工具箱，用 MATLAB勺FDATool工具设计该滤波器，参数设置如图5所示：设计好参数后，点击DesignFilter，可以在FDATool窗口的左上角看到滤波器的阶数为63阶，点击 File Gen erate M-file ，并将滤波器命名为 mylowfilter 。D y # 妙図 D E田誥击阳J- M O K 0 皙l Current Fitter rilormation一 MetgnrtiJcte Response (dBjStructurn:Diredt-Formi FIROrder:63StbUfeYesSoMrce：DesignedStore

9、Filer.Fiter lilanger(Response TpeHtghpassBafbdpassBafidstapDifferentiator Design MethodIR |QutterwgTth11-800 2Freqiiemcy (MHz10 12低通滤波器(Filler OrderC Specify orderf* I inirrum zrderl OptiansDensity Factor. |H金 Fin jEquiripple寸0频率采样中的等兹纹舞壘Unts: dBApess：A slop；采样频25MHzl Magnitude Specllications.BO阻带袞减

10、SOdEt|Desigrting Filler . Done图5编写如下代码，滤除混叠信号中的高频信号。Hd = mylowfilter; %滤波器名称output = filter(Hd, y1); %对叠加信号y，进行滤波处理figure(3);plot(y2, k);hold on;plot(output, r);legend(0.5MHz 原始信号,滤波后取出的信号);title(信号通过MTALAB的低通滤波器后的波形);滤波后的波形如图6所示:生成.coe文件，用于Xilinx 的IPCore设计滤波器，FDATool窗口点击File Export, 保持默认设置，点击 Expor

11、t即可，次数在 MATLAB勺workspace窗口多出一个 Num的1*64 的数组，这就是滤波器的系数，如图7所示：fals d】t Aulji i a Ticgh I i Vi hv liin Hip宙QQX O E S k O 耕创西钊足 slut K三| 3*11(2*1lBFittei CHiiidenifD.-K4z24EICiX9i4 5：-!-!2E4ei D.OO4-HG4S；a FTX3T-S5I3Z D .47a177&Ca-EE-M7J931!1 Q W 中7理朝l曰测D.-Mi23fii3ViieMJ43?.0&133BB7B5!I32D.COL13SK77D4E.

12、M id-&破射加初】*治！H科甲门-.-M-i2flifi3aT!：5Ma?231HalpcfiAi IVP*. FCnllHtaillQIF-l_=A；AiSXlaxBri IMhcd.r 詔 iHunsTA-crti洽 5 |Eqnip4eUnite |UH1二LHt 肚jJFa：旷 伽 5： |1*fM- fAzlafr 脚FSKfr k_Fvi Lrdn-F Spacey 3 nferCompirlinq RcsporiM ctone打品H肌*1305WKrWWW3BLBSOW2WKSWOGL他召ClrfKJ AnHa?书 eflLXKOUxfilXXl dnila-fiOSKrJ

13、MJ.a-D ：aO.M4.L&04.gClrfKJ dnbla333057srsu tbt ISlaI I *33767) % klni：i 玮:期td，. fii b fqpwi r #； /f irrpr t?i+r i!:由于MATLABfc成的滤波器系数全是一些小数，而FPGA只能处理整数，因此我们必须将这些小数扩大一定的倍数，使它们变成整数。在MATLAB勺命令窗口输入下面的一段代码然后按Enter，即可将上面这些系数变为整数，注意这里的*32767，表示将系数扩大 32767倍，这里的扩大倍数只能选 2AN,目的是为 了后面滤波后的波形数据的高位截取（丢掉低位，即除以2AN）。返

14、回ans=0,表示操作正确。coeff=rou nd(Num/max(abs(Num)*32767);% abs() 求绝对值，max()求最大值，round() 四舍五入fid = fope n(e:/fircoe.txt,wt); %将滤波器系数写入文件件中fprintf(fid,%16.0fn,coeff);%将滤波器系数以16位浮点数的格式保存fclose(fid)程序运行的结果如图8所示：将文件的格式改为.coe格式，在文件的开口加上:radix = 10;coefdata =在每个数字前的空格处加一个逗号“，”，并在最后一个数字后面加上分号“；”如图9所示：12r地1Jr14 ,小

15、 15 (-1544起as住诚号TITi12r-4EZ&1?r2&-树姑.-S12,- = ；: =23 匚-151724 |-2515-这样就产生了可供 FIR滤波器的IPCore可调用的.coe文件。3）用Verilog编程读取MATLAB生的波形数据通过Verilog 的$readmemh（）函数即可读出.txt 文件中的数据（这里的数据是以16进制形式存储的），具体用法如下：reg 15 :0 data_i n0:2000; /定义一个 16bit * 2001 的数组initial begin读出 MATLAB生的波形数据 0.5MHz_sin + 4MHz_sin 信号$readm

16、emh（c:/matlab_wave_data.txt, data_ in）;/将matlab_wave_data.txt中的数据读入存储器data_inEnd我在modelsim中尝试用readmemh来读取matlab生成的数据，但是一旦到了负数的地方就出现问题了，之后的所有数据都变成了xxxxx，因此我在用 MATLAB生波形数据时都给正弦信号添加了一个偏移量，将正弦信号平移到了0之上。用Modelsim仿真工具可以看到读出的数据，并与原文件数据（图 3）比较，可以判断读出的数据是否正确，这里是ISE和Modelsim联合仿真的结果，如图 10所示：图104） 由Xilinx 的FIR的

17、IPCore实现FIR滤波器FIR的IPCore的生成步骤就不多说了，简单的说一下参数配置情况：a）我们使用由MATLABfc成的.coe文件导入滤波器的系数。b）输入信号的采样频率为25MHz这里必须和MATLAB产生的波形数据的采样频率保持一致。c）FIR滤波器的运行时钟，这里我们选择250MHz这里必须和仿真文件里给的时钟保持一致。d）输入数据的类型，我们这里选择了无符号数，输入数据位宽为16bit，因为modelsim不能读出负数（要么用补码，要么加偏移）。如图11、图12所示：，-I Ml.rfd(rfd), / output rfd在其上升沿将输入数据加载到滤波器内核中Freq.

18、ResponseFIR CompilerlrM.ciKTi:ip:lir_mpilr75LOFrequency Response (Magnitude)iHDfTFialLZd F僂4血吟(k FQ rad/豈汙电1：FJisr Sf*ew$n輸.扎舊号采徉频率25MII上QrtBShMtenenrie |Cancel |Helpidtk*秤 1 of-*邛lay : iRnngc: 1.1Aber ArhsFpgzdndSbop bandHant: 10.0- iLS |0.S- | L.OMn :-2&.Z52fi44 cE融H :LiO】M昭曲35.W27+援Ripe : LID.UTS

19、DUdBMPSyflibd .所*牛和呷5骑| Q ImpkmLMcm血图11图12FIR的IPCore的列化如下：FIR16_IP FIR16_IP_i ns (.clk(clk), /in put elk.rdy(rdy), / output rdy在其上升沿输出滤波器的计算结果=0;end特别注意这个数据位宽我们主要对其进行简单的控制:在rfd上升沿将输入数据加载到滤波器内核中,在rdy上升沿输出滤波器的计算结果。具体的Verilog 代码如下:.di n(data_in_reg), / i nput 15 : 0 din.dout(dout); / output 35 : 0 dout

20、always (posedge elk) begi n信号的上升沿将输入数据加载到滤波器内核中if(reset = 1b0) begi ni=1b0;m=1b0;data_in_reg=16h0000;beginrfd_1q =rfd;rfd_2q =rfd_1q;endif(rfd_1q & rfd) begin / rfdelseendendenddata_i n_reg=if(i = 2002)always (posedge elk) begi nif(reset = 1b0) begi nData_out_reg=data_ in i;i + 1;m;0;0;0;=else begi

21、nrdy_1q =rdy;rdy_2q =rdy_1q;if(rdy_1q & rdy) begi nendendendData_out_reg=dout;j + 1b1;n;这里还做了一个附加功能，将FIR滤波器的输入数据存放到一个.txt文件当中，然后用MATLAB读取这个波形文件数据,看看读出的波形是否和原来的混叠波形一样。具体的Verilog和MATLAB弋码如下:in tegerwr_file;in itialwr_file = $fope n(c:/FIR_in_data.txt);always ( m ) begi nif(reset = 1b1) begi n$fdisplay

22、(wr_file,%h, data_in_reg);/ 33bitif(j = 11d2002)/共写入 2001 个数据$stop;end end fid = fope n( c:/FIR_in_data.txt,丫);这句话的意思是从fid所指的文件以16进制方式读出for j = 1 : 2000;num1(j) = fscan f(fid, %x, 1);%一个数据。end fclose(fid);figure(4);plot( nu ml,丫)；legend(Verilog 读出的txt文件中的数据);title(FIR 滤波器的输入数据);MATLAB读出的波形数据如图 13所示:

23、图135)将FIR滤波器的输出存入.txt文件中有上面的步骤我们已经验证了输入到FIR滤波器的数据都是正确的，下面就将FIR滤波器的输出数据保存到一个.txt文件当中供MATLAB卖取。/刚刚的问题是， matlab读一个数据是32bit的，/而FIR的输出是36bit的，因此高4bit根本没有读上来。assig nData_out31:0 = Data_out_reg35:4;in tegerw_file;in itialw_file = $fope n(c:/FIR_out.txt);always ( n ) beginif(reset = 1b1) begi n$fdisplay(w_f

24、ile,%h, Data_out31:14); 33bit数if(j = 11d2002)/共写入 2001 个数据$stop;endend这里我也搞了好久才搞好，这里FIR滤波器的输出数据位宽变成了36bit ,而输入数据位宽是16bit，为什么数据会变大几万倍呢？因为我们在将滤波器的系数由小数变成整数的时候，对这些系数整体扩大了32767倍,再做了一个四舍五入（影响滤波器精度），对滤波器的系数扩大的倍数越大，四舍五入对精度的影响就越小，但是系数乘的倍数越大，FPGA在做乘加运算也就越复杂，也就越耗时，越耗资源，因此我们需要找一个平衡点。这里为了将信号的幅度变回原始的幅度（尽可能的靠近），我

25、们只能通过将低位截取掉，截取低位相当于对数据做除法（除2），所以前面的滤波器系数的扩大倍数我们一定要用2AN，这样我们在这里还原信号幅度的时候，只需要截位就能达到目的。比如这里我们对滤波器的系数乘了32767，那我们在做除法还原波形幅度时，只需要除以32767即可（即截掉低16bit ）。还有一种操作方式就是我们只保留数据的高16bit （和输入数据的位宽保持一致），这两种方式波形的幅度也就几倍的差距，我还没有完全搞懂这里，究竟怎样才能将波形的幅度完全的还原回去，还是一个值得好好思考的问题？另一个问题是我的电脑是 32bit位宽的，如果我们一次性让 MATLAB取36bit的数据 那么数据的高

26、4bit会读不上来，会导致很奇怪的波形，我也遇到了这个问题。 如图14所示， 波形明显可以通过一些平移拼合成一个正弦波。通过对波形数据一个个的分析， 我找到了这个问题。图14图15是波形幅度发生变化的截图:图156）由MATLAB入FIR滤波器的输出数据，并分析滤波结果由MATLAB取FIR滤波器的输出数据，分析波形，具体的代码如下:fid = fope n( c:/FIR_out.txt,丫);for i = 1 : 2000;num(i) = fscanf(fid, %x, 1);%这句话的意思是从fid所指的文件以16进制方式读出一个数据。endfclose(fid);figure(5)

27、;plot (n um,丫)；y4 = y2;hold on;plot(y4, k);legend(经过FIR_IPCore滤波后的数据,0.5MHz的原始数据放大16384倍);title( 经过FIR滤波器的输出数据);波形数据如图16所示:图16总结本次设计我大约耗时一周，从0开始研究FIR滤波器的设计，联合了 MATLAB ISE、Modelsim三个工具，不得不佩服 MABLA莊数据分析方面强大的功能。本次设计走通了 FIR滤波器的总体设计流程，为以后的工程实用打下了基础，当然还有一些基本问题， 如波形幅度如何完美的还原还没搞清楚，有待进一步研究。在实际工程应用时，我们可以通过联合

28、MATLAB ISE、Modelsim 三个工具完全脱离硬件来设计、仿真、验证FIR数字滤波器的性能，肯定能够大大的缩短设计周期，提升滤波器性能。附录：附录为Verilog源代码和MATLAB源代码，这些源代码是经过调试的，是可以直接使用的。供大家参考Verilog 源代码：module FIR_Lowpass( clk, reset, Data_out);in putclk;in putreset;output31 :0 Data_out;reg35 :0 Data_out_reg;reg 10:0 i = 0;reg 10:0 j = 0;reg 15 :0 data_i n0:2000;

29、 /定义一个 16bit * 2001 的数组reg 15:0 data_ in_reg = 0;initial begin读出 MATLAB生的波形数据 0.5MHz_sin + 4MHz_sin 信号$readmemh(c:/matlab_wave_data.txt, data_i n); /将 matlab_wave_data.txt 中的数据读入存储器data_i nendwirerfd;wirerdy;wire 35:0 dout;reg rfd_1q;reg rfd_2q;reg rdy_1q;reg rdy_2q;reg n = 0;reg m = 0;always (posed

30、ge clk) begi nif(reset = 1b0) begi ni=1b0;m=1b0;data_in_reg=16h0000;endelse begi nrfd_1q=rfd;rfd_2q=rfd_1q;if(rfd_1q & rfd) begin / rfd信号的上升沿将输入数据加载到滤波器内核中data_in_regi m=if(i = 2002)iendendendFIR16_IP FIR16_IP_i ns (.clk(clk), /in put elk .rfd(rfd), / output rfd .rdy(rdy), / output rdy=data_i ni;=i

31、+ 1;m;=0;在其上升沿将输入数据加载到滤波器内核中在其上升沿输出滤波器的计算结果.di n(data_in_reg), / i nput15 : 0 din.dout(dout); / output 35 : 0 doutalways (posedge elk) begi nif(reset = 1b0) begi nData_out_reg=0;0;=0;=endelse begi nrdy_1q =rdy;rdy_2q =rdy_1q;if(rdy_1q & rdy) begi nendendendData_out_reg=dout;j + 1b1;n;/刚刚的问题是，matlab读

32、一个数据是 32bit的，而FIR的输出是36bit的，因此高4bit根本没有读上来。/ assig n Data_out31:0 = Data_out_reg35:4;assig nData_out31:0 = Data_out_reg35:4;/*in tegerwr_file;in itialwr_file = $fope n(c:/FIR_in_data.txt);always ( m ) begi nif(reset = 1b1) begi n$fdisplay(wr_file,%h, data_in_reg);/ 33bit数if(j = 11d2002)/共写入 2001 个数据

33、$stop;endend/*in tegerw_file;in itialw_file = $fope n(c:/FIR_out.txt);always ( n ) beginif(reset = 1b1) begi n$fdisplay(w_file,%h, Data_out31:14); 33bit数if(j = 11d2002)/共写入 2001 个数据$stop;endend en dmoduleMATLAB源代码/6*产生信号并保存到.txt文件中*clear all;fs = 25000000; % 25M 采样率t = 0:1/fs:0.0002;%共 0.0002 * 2500

34、0000 = 5000 个点fl = 500000;f2 = 2000000;signal1 = sin(2*pi*f1*t); %频率为0.5MHz的正弦信号sig nal2 = sin( 2*pi*f2*t); %频率为4.0MHz的正弦信号%y1 = signal1 + signal2; %两个正弦信号叠加%x = lin space(0, 12.56, 2048); %在区间0,6.28 = 2*pi之间等间隔的取1024个点%y1 = sin (x); %计算相应的余弦值%由于正、余弦波形的值在0,1之间，需要量化成16bit,先将数值放大 %y1 = y1 * 32768; % 3

35、2 * 1024 = 32768 %y1 = y1 * 16384; % 32 * 1024 = 32768%y1 = y1 + 32768;y2 = fix(16384 + (2A14 - 1) * sig nal1);y3 = fix(16384 + (2人14 - 1) * signal2);y1 = y2 + y3;%再将放大的浮点值量化，并写道到存放在C盘的文本中fid = fope n(c:/matlab_wave_data.txt, wt);%fprintf(fid, %16.0fn, y1); %在写文件时量化为16bit 的定点实数【%16.0f,16.0表示16bit定点数，f表示实数】，范围是：-32768 - 32767fprintf(fid, %xn, y1); %在写文件时量化为16bit 的定点实数【%16.0f,16.0 表示16bit定点数，f表示实数】，范围是：-32768 - 32767fclose(fid);figure(1);plot(y2, b);hold on;plot(y3, b);hold on;plot(y1, r);legend（0.5MHz 正弦,2MHz 正弦,两者叠加）; title（MATLAB 产生的两个正弦信