基于FPGA的16階級(jí)聯(lián)型iir帶通濾波器實(shí)現(xiàn)

goandlove 2019-09-11

展開全文

基礎(chǔ)知識(shí)：

什么叫濾波器？

簡單的說，就像篩米，留下你需要的米，濾掉不需要的米頭。過濾的功能。

什么叫數(shù)字濾波器？

用數(shù)字芯片做的濾波器，而不是rc搭的，輸入是離散的序列，輸出也是離散的序列；

快速了解時(shí)域頻域：

https://zhuanlan.zhihu.com/p/19763358?from=singlemessage&isappinstalled=1

什么叫時(shí)域？

信號(hào)隨時(shí)間的變化。

什么叫頻域？

曾經(jīng)有個(gè)通俗的解釋是：彈鋼琴，琴鍵1234等表示的就是頻域，產(chǎn)生的各種音樂就是時(shí)域，你以為的萬變其實(shí)是永恒的不變。

什么叫fir與iir濾波器？

FIR(Finite Impulse Response)濾波器：有限長單位沖激響應(yīng)濾波器，又稱為非遞歸型濾波器，是數(shù)字信號(hào)處理系統(tǒng)中最基本的元件，它可以在保證任意幅頻特性的同時(shí)具有嚴(yán)格的線性相頻特性，同時(shí)其單位抽樣響應(yīng)是有限長的，因而濾波器是穩(wěn)定的系統(tǒng)。

無限脈沖響應(yīng)。遞歸濾波器，也就是IIR數(shù)字濾波器，顧名思義，具有反饋。

fir和iir有啥異同（important）？

根據(jù)沖激響應(yīng)的不同，將數(shù)字濾波器分為有限沖激響應(yīng)（FIR）濾波器和無限沖激響應(yīng)（IIR）濾波器。對(duì)于FIR濾波器，沖激響應(yīng)在有限時(shí)間內(nèi)衰減為零，其輸出僅取決于當(dāng)前和過去的輸入信號(hào)值。對(duì)于IIR濾波器，沖激響應(yīng)理論上應(yīng)會(huì)無限持續(xù)，其輸出不僅取決于當(dāng)前和過去的輸入信號(hào)值，也取決于過去的信號(hào)輸出值。

1. 在相同技術(shù)指標(biāo)下，IIR濾波器由于存在著輸出對(duì)輸入的反饋，因而可用比FIR濾波器較少的階數(shù)來滿足指標(biāo)的要求，這樣一來所用的存儲(chǔ)單元少，運(yùn)算次數(shù)少，較為經(jīng)濟(jì)。例如用頻率抽樣法設(shè)計(jì)阻帶衰減為-20db的FIR濾波器，其階數(shù)要33階才能達(dá)到，而如果用雙線性變換法設(shè)計(jì)只需4-5階的切貝雪夫?yàn)V波器，即可達(dá)到指標(biāo)要求，所以FIR濾波器的階數(shù)要高5-10倍左右。

2. FIR濾波器可得到嚴(yán)格的線性相位，而IIR濾波器則做不到這一點(diǎn)，IIR濾波器選擇性愈好，則相位的非線性愈嚴(yán)重，困而，如果IIR濾波器要得到線性相位，又要滿足幅度濾波的技術(shù)要求，必須加全通網(wǎng)絡(luò)進(jìn)行相位校正，這同樣會(huì)大大增加濾波器的階數(shù)，從這一點(diǎn)上看，F(xiàn)IR濾波器又優(yōu)于IIR濾波器。

3. FIR濾波器主要采用非遞歸結(jié)構(gòu)，因而從理論上到實(shí)際的有限精度的運(yùn)算中，都是穩(wěn)定的。有限精度運(yùn)算誤差也較小，IIR濾波器必須采用遞歸的結(jié)構(gòu)，極點(diǎn)必須在Z平面單位圓內(nèi)，才能穩(wěn)定，這種結(jié)構(gòu)，運(yùn)算中的四舍五入處理，有時(shí)會(huì)引起寄生振蕩。

4. FIR濾波器，由于沖激響應(yīng)是有限長的，因而可以用快速傅里葉變換算法，這樣運(yùn)算速度可以快得多，IIR濾波器則不能這樣運(yùn)算。

5. 從設(shè)計(jì)上看，IIR濾波器可以利用模擬濾波器設(shè)計(jì)的現(xiàn)成閉合公式、數(shù)據(jù)和表格，因而計(jì)算工作量較小，對(duì)計(jì)算工具要求不高。FIR濾波器則一般沒有現(xiàn)成的設(shè)計(jì)公式，窗函數(shù)法只給出窗函數(shù)的計(jì)算工式，但計(jì)算通帶、阻帶衰衰減仍無顯示表達(dá)式。一般FIR濾波器設(shè)計(jì)只有計(jì)算機(jī)程序可資利用，因而要借助于計(jì)算機(jī)。

6. IIR濾波器主要是設(shè)計(jì)規(guī)格化的、頻率特性為分段常數(shù)的標(biāo)準(zhǔn)低通、高通、帶通、帶阻、全通濾波器，而FIR濾波器則要靈活得多，例如頻率抽樣設(shè)計(jì)法，可適應(yīng)各種幅度特性的要求，因而FIR濾波器則要靈活得多，例如頻率器可設(shè)計(jì)出理想正交變換器、理想微分器、線性調(diào)頻器等各種網(wǎng)絡(luò)，適應(yīng)性較廣。而且，目前已有許多FIR濾波器的計(jì)算機(jī)程序可供使用。

什么叫定點(diǎn)數(shù)？

計(jì)算機(jī)中采用的一種數(shù)的表示方法。參與運(yùn)算的數(shù)的小數(shù)點(diǎn)位置固定不變。

什么叫濾波器的零點(diǎn)極點(diǎn)？

濾波器可以看成是一個(gè)信號(hào)處理的系統(tǒng),其輸入輸出之間存在一定的關(guān)系,這種關(guān)系無論在時(shí)域還是頻域都可以用數(shù)學(xué)表達(dá)式來表示.而這數(shù)學(xué)表達(dá)式又是分子分母都是多項(xiàng)式的表達(dá)式（稱為傳輸函數(shù)）,這樣滿足使傳輸函數(shù)的分子為零的是零點(diǎn),滿足使傳輸函數(shù)分母為零的就是其極點(diǎn).

iir濾波器的種類：很多啊，直接一型，直接二型，級(jí)聯(lián)型，并聯(lián)型。

對(duì)于matlab的fdatool工具中二階節(jié)默認(rèn)結(jié)構(gòu)為：

對(duì)于這個(gè)結(jié)構(gòu)用圖表示為：

差分方程表示為：

零極點(diǎn)表示為：零點(diǎn)就是差分方程的前面三項(xiàng)，極點(diǎn)就是后面兩項(xiàng)。用FPGA實(shí)現(xiàn)主要就是實(shí)現(xiàn)濾波器的差分方程。

流程：

任務(wù)要求：

16階二階級(jí)聯(lián)IIR數(shù)字濾波器設(shè)計(jì)，16bit有符號(hào)整數(shù)連續(xù)輸入，采樣率80khz，通帶頻率1k-8khz。系數(shù)為16bit有符號(hào)整數(shù)。

1.系數(shù)產(chǎn)生：通過matlab中的fdatool軟件生成所需系數(shù)。（當(dāng)然可以用各種函數(shù)生成，太難工科生表示要陣亡了，還是默默用fdatooll吧）

把需求放入fdatool中：生成的架構(gòu)就是直接二型二階節(jié)結(jié)構(gòu)。

零極點(diǎn)圖：

未量化的系數(shù)：

未量化的系數(shù)導(dǎo)出：生成一個(gè)c文件。

那么問題來了，這個(gè)c文件中的內(nèi)容是啥子意思呢，一開始我也是一臉懵逼，而且網(wǎng)上的資料少之又少，文件如下所示，含義已注釋：

1 /*
2 * Filter Coefficients (C Source) generated by the Filter Design and Analysis Tool
3 *
4 * Generated by MATLAB(R) 7.8 and the Signal Processing Toolbox 6.11.
5 *
6 * Generated on: 22-Sep-2017 20:23:35
7 *
8 */
9
10 /*
11 * Discrete-Time IIR Filter (real)
12 * -------------------------------
13 * Filter Structure : Direct-Form II, Second-Order Sections
14 * Number of Sections : 8
15 * Stable : Yes
16 * Linear Phase : No
17 */
18
19 /* General type conversion for MATLAB generated C-code */
20 #include 'tmwtypes.h'
21 /*
22 * Expected path to tmwtypes.h
23 * D:\workfile\Matlab2009\extern\include\tmwtypes.h
24 */
25 /*
26 * Warning - Filter coefficients were truncated to fit specified data type.
27 * The resulting response may not match generated theoretical response.
28 * Use the Filter Design & Analysis Tool to design accurate
29 * single-precision filter coefficients.
30 */
31 #define MWSPT_NSEC 17
32 const int NL[MWSPT_NSEC] = { 1,3,1,3,1,3,1,3,1,3,1,3,1,3,1,3,1 };
33 //上面1313的玩意表示下面這個(gè)數(shù)組哪個(gè)項(xiàng)有效，1則表示第一項(xiàng)有效，3表示都有效；
34 const real32_T NUM[MWSPT_NSEC][3] = {
35 {
36 0.1001105756, 0, 0 //第一個(gè)二階節(jié)的增益；
37 },
38 {
39 1, 0.7806397676, 1 //第一個(gè)二階節(jié)的零點(diǎn)；b0,b1,b2;
40 },
41 {
42 0.1001105756, 0, 0 //第二個(gè)二階節(jié)的增益；
43 },
44 {
45 1, -1.999714136, 1 //第二個(gè)二階節(jié)的零點(diǎn)；b0,b1,b2;
46 },
47 {
48 0.3725369573, 0, 0 //以下就是類似的了；
49 },
50 {
51 1, -0.9795594215, 1
52 },
53 {
54 0.3725369573, 0, 0
55 },
56 {
57 1, -1.99809742, 1
58 },
59 {
60 0.6452683806, 0, 0
61 },
62 {
63 1, -1.352879047, 1
64 },
65 {
66 0.6452683806, 0, 0
67 },
68 {
69 1, -1.996625185, 1
70 },
71 {
72 0.7896357179, 0, 0
73 },
74 {
75 1, -1.448690891, 1
76 },
77 {
78 0.7896357179, 0, 0
79 },
80 {
81 1, -1.995926261, 1
82 },
83 {
84 1, 0, 0 //總的增益為1，上面8個(gè)分增益相乘最終為1；
85 }
86 };
87 const int DL[MWSPT_NSEC] = { 1,3,1,3,1,3,1,3,1,3,1,3,1,3,1,3,1 };
88 const real32_T DEN[MWSPT_NSEC][3] = {
89 {
90 1, 0, 0 //忽略項(xiàng)；
91 },
92 {
93 1, -1.765431523, 0.8537048697 //第一個(gè)二階節(jié)的極點(diǎn)；a0,a1,a2;
94 },
95 {
96 1, 0, 0
97 },
98 {
99 1, -1.893844962, 0.919323802 //以下類似；
100 },
101 {
102 1, 0, 0
103 },
104 {
105 1, -1.666594863, 0.877212882
106 },
107 {
108 1, 0, 0
109 },
110 {
111 1, -1.959967136, 0.9707458019
112 },
113 {
114 1, 0, 0
115 },
116 {
117 1, -1.614711642, 0.9346644878
118 },
119 {
120 1, 0, 0
121 },
122 {
123 1, -1.982463837, 0.9896451831
124 },
125 {
126 1, 0, 0
127 },
128 {
129 1, -1.603200555, 0.9806866646
130 },
131 {
132 1, 0, 0
133 },
134 {
135 1, -1.991223216, 0.9973948002
136 },
137 {
138 1, 0, 0
139 }
140 };

系數(shù)量化選項(xiàng)：系數(shù)量化你可以自己量化也可以讓軟件量化，不過它量化出來的數(shù)據(jù)零點(diǎn)并不是乘完增益后再進(jìn)行量化的。最好還是乘完增益后再量化，所以還是自己用excel慢慢量化吧，眼淚掉下來。

未量化excel表：

excel中計(jì)算單元格方便到不行：零點(diǎn)乘完增益放大16384；極點(diǎn)直接放大16384；下圖gain請(qǐng)無視。

新的b0=b0*gain1*16384;新的a0=a0*16384;放大16384倍方便FPGA實(shí)現(xiàn)除法截位。

2.編碼實(shí)現(xiàn)：

先看一下16階iir濾波器架構(gòu)：級(jí)聯(lián)8個(gè)二階節(jié)。

一個(gè)二階節(jié)：

現(xiàn)在就可以編碼實(shí)現(xiàn)它了，這是第一版代碼，尚未優(yōu)化，仿真ok，不要邏輯綜合，會(huì)占用成噸的資源。

由于技術(shù)垃圾，不做十分精確輸出位控制，輸出都為16bit數(shù)據(jù)。

兩個(gè)n位的加法結(jié)果需要n+1位；兩個(gè)n位的乘法結(jié)果需要2n位。

matalb生成modelsim仿真文件向量：

生成1500hz，采樣80khz波形向量文件。生成其他hz的波形文件類似。

 1 f1=1500;   %頻率1500hz；
 2 Fs=80000;  %采樣80khz；
 3 N=16;        %16bit量化；
 4 t=0:1/Fs:0.01;  %采樣時(shí)長0.01；
 5 c2=2*pi*f1*t;
 6 s2=sin(c2);  %正弦波產(chǎn)生；
 7 s2=s2/max(abs(s2));
 8 Q_s=round(s2*(2^(N-1)-1));
 9 plot(t,s2,'r*-');   %畫圖；
10 
11 fid=fopen('D:\data\data_1500\data_1500.txt','w');    %采樣點(diǎn)保存為10進(jìn)制；
12 fprintf(fid,'%8d\r\n',s2);
13 fprintf(fid,';');
14 fclose(fid);
15 
16 fid=fopen('D:\data\data_1500\data_1500_B.txt','w'); %采樣點(diǎn)保存為2進(jìn)制；
17 for i=1:length(Q_s)
18     B_s=dec2bin(Q_s(i)+(Q_s(i)<0)*2^N,N)
19     for j=1:N
20        if B_s(j)=='1'
21            tb=1;
22        else
23            tb=0;
24        end
25        fprintf(fid,'%d',tb);
26     end
27     fprintf(fid,'\r\n');
28 end
29 fprintf(fid,';');
30 fclose(fid);

仿真測試：

對(duì)600hz正弦波濾波結(jié)果：600hz波形被濾除。

對(duì)5000hz正弦波濾波結(jié)果：5000hz波形通過。

對(duì)9000hz波形濾波結(jié)果：開始有點(diǎn)點(diǎn)迷之振蕩，基本濾除9000hz的波。

最開始的結(jié)果經(jīng)過多久出來到out？（特么上次面試還問這個(gè)了，十臉懵逼，根本沒注意這啊。。。emmm很氣）

可以看到是復(fù)位拉高后的9個(gè)時(shí)鐘周期后yout數(shù)據(jù)產(chǎn)生，因?yàn)榱魉€啊，emmm。

初版代碼綜合上板子：通過rom輸出5khz的數(shù)據(jù)。

所以優(yōu)化很重要，這是未優(yōu)化版本。

signaltapII抓下波：

優(yōu)化版以及未優(yōu)化版比較：只包含iir部分，不含pll以及rom。系統(tǒng)時(shí)鐘跟采樣時(shí)鐘一樣，80khz。

未優(yōu)化版：直接采用*（乘）的方式。

優(yōu)化版：采用內(nèi)置乘法器，以及采用移位相加的方法。資源少的可憐啊，一共才30個(gè)9bit乘法器。。。。，若再增加乘法器，le使用量又會(huì)往上漲。未來優(yōu)化方向：提高時(shí)鐘頻率，復(fù)用乘法器。

其他：

怎么優(yōu)雅的分解系數(shù)用來移位相加：

直接寫了個(gè)c程序，來看看效果：

c源代碼：看看就好啦，很久沒寫c，完全沒有代碼style了emmm。

1 #include <stdio.h>
2 #include <math.h>
3 int main(void)
4 {
5 int coefficient;
6 int sum;
7 int sum1;
8 int mul;
9 int mul1;
10 int j;
11 int i;
12 int k=0;
13 int m;
14 int n=0;
15 int cha;
16 printf('All rights by kingstacker!\n');
17 begin:
18 printf('Pelese input the coefficient:');
19 scanf('%d',&coefficient);
20 printf('%d=',coefficient);
21 sum = coefficient;
22 sum1 = coefficient;
23 for (m=15;m>=0;m--) //add;
24 {
25 mul1=pow(2,m);
26 if (sum1 >= mul1)
27 {
28 sum1 = sum1 -mul1;
29 n=n+1;
30 printf('+%d(2^%d)',mul1,m );
31
32 }
33
34 }
35 printf('\nIf add,use %d add source !\n',n-1 );
36 //sub;
37 for (j=0;j<=15;j++)
38 {
39 mul=pow(2,j);
40 if (mul >= sum)
41 {
42 goto this;
43 }
44 }
45 this:
46 cha = mul - sum;
47 printf('%d=%d(2^%d)',sum,mul,j );
48 for (i=j;i>=0;i--)
49 {
50 mul1 = pow(2,i);
51 if (cha >= mul1)
52 {
53 cha = cha - mul1;
54 k=k+1;
55 printf('-%d(2^%d)',mul1,i );
56 }
57 }
58 printf('\nIf sub,use %d add source !\n',k );
59 //result;
60 if((n-1) <= k)
61 {
62 printf('\nadd is better!\n');
63 }
64 else
65 {
66 printf('\nsub is better!\n');
67 }
68 k=0;
69 n=0;
70 goto begin;
71 printf('Thanks for you use!bye!\n');
72
73 }

以上。

歡迎通信工程師和FPGA工程師關(guān)注公眾號(hào)

FPGA微信技術(shù)群

歡迎大家加入全國FPGA微信技術(shù)群，這里有一群熱愛技術(shù)的工程師，在這里可以一起交流討論技術(shù)！

用手指按住就可以加入FPGA全國技術(shù)群哦

FPGA IP核服務(wù)：各類優(yōu)質(zhì)IP核服務(wù)商，服務(wù)到位，有保障！有需求的可以直接聯(lián)系群主！

本站是提供個(gè)人知識(shí)管理的網(wǎng)絡(luò)存儲(chǔ)空間，所有內(nèi)容均由用戶發(fā)布，不代表本站觀點(diǎn)。請(qǐng)注意甄別內(nèi)容中的聯(lián)系方式、誘導(dǎo)購買等信息，謹(jǐn)防詐騙。如發(fā)現(xiàn)有害或侵權(quán)內(nèi)容，請(qǐng)點(diǎn)擊一鍵舉報(bào)。

轉(zhuǎn)藏 分享

QQ空間 QQ好友新浪微博微信

獻(xiàn)花（0） +1

來自： goandlove > 《信號(hào)處理》

舉報(bào)/認(rèn)領(lǐng)