基于FPGA IIR數(shù)字濾波器的設計

IIR（Infinite Impulse Response）無線脈沖響應濾波器。
系統(tǒng)傳遞函數(shù)為：

系統(tǒng)的差分方程可寫為：

IIR優(yōu)缺點:
1)在相同的幅頻條件下，濾波器階數(shù)比FIR濾波器低。
2)IIR濾波器占用的硬件資源比較少（相比FIR濾波器）。
3)不具備嚴格的線性相位特性。

1 IIR數(shù)字濾波器的基本結構及類型

圖1 直接I型

圖2 直接II型

圖3 級聯(lián)型

圖4 并聯(lián)型

2 設計目標

采用matlab buffer函數(shù)設計一個IIR濾波器低通濾波器，通帶截止頻率為1khz，輸入信號為1khz 3khz sin波形，經過IIR濾波器后輸出為1KHZ sin波，其他不做要求。（本文只對IIR設計思想進行驗證不做性能要求）。

3 matlab的設計驗證

Matlab源碼：

%參數(shù)定義

FS =44100; %Sample rate Frequncy
fc = 1000; %1khz
fe = 3000; %外部輸入信號 3khz
N = 1024;
Q =16;
%波形產生

sin_osc =sin(t*fc);
sin_e =sin(t*fe);
sin_add = sin_osc sin_e;
%IIR 濾波器系數(shù)（低通濾波器）

[b a] = butter(3,fc/(FS/2),'low');

%濾波（混頻后）
y = filter(b,a,sin_add);

f_osc =fft(sin_osc,N);
f_osc=20*log(abs(f_osc))/log(10); %換算成dBW單位
ft=[0:(FS/N):FS/2]; %轉換橫坐標以Hz為單位
f_osc=f_osc(1:length(ft));
%濾波器系數(shù)量化

Mab =max(max(abs(a),abs(b)));
%16bit 量化
Qb = round((b/Mab)*(2^(Q-1)-1));
Qa = round((a/Mab)*(2^(Q-1)-1));
%%%%
Qm =floor(log2(Mab/a(1)));
if Qm<log2(Mab/a(1))
Qm = Qm 1;
end
Qm = 2^Qm;
Qb1=round(b/Qm*(2^(Q-1)-1));
Qa1=round(a/Qm*(2^(Q-1)-1));
%繪圖
%時域波形圖

figure(1),
hold on
subplot(221),plot(t(1:128),sin_osc(1:128),'-');
legend('sin 1khz');title('sin 1khz');
subplot(222),plot(t(1:128),sin_e(1:128),'-');
legend('sin 3khz');title('sin 3khz');
subplot(223),plot(t(1:128),sin_add(1:128),'-');
legend('sin 1khz add 3khz');title('sin 1khz add 3khz');
subplot(224),plot(t(1:128),y(1:128),'-');
legend('LPF 結果');title('LPF 結果');
grid;
hold off
%頻域波形

figure(2),
hold on
subplot(221);plot(ft,f_osc);
xlabel('頻率(Hz)','fontsize',8); ylabel('功率(dBW)','fontsize',8);
title('信號頻譜圖 2KHZ','fontsize',8);legend('sinosc');
subplot(222);plot(ft,f_e);
xlabel('頻率(Hz)','fontsize',8); ylabel('功率(dBW)','fontsize',8);
title('信號頻譜圖3KHZ','fontsize',8);legend('sine');
subplot(223);plot(ft,f_add);
xlabel('頻率(Hz)','fontsize',8); ylabel('功率(dBW)','fontsize',8);
title('信號頻譜圖2KHZ 和 3KHZ','fontsize',8);legend('sin add');
subplot(224);plot(ft,y_f);
xlabel('頻率(Hz)','fontsize',8); ylabel('功率(dBW)','fontsize',8);
title('信號頻譜圖濾波后','fontsize',8);legend('LPF結果');
hold off
%幅頻響應

figure(3);
subplot(211);stem(Fb);
title('Fb單位抽樣響應','fontsize',8);
subplot(212);plot(f,mag);
xlabel('頻率(Hz)','fontsize',8);
ylabel('幅度(dB)','fontsize',8);
title('freqz()幅頻響應','fontsize',8);

圖5 matlab時域波形如圖5所示，（5,1）波形加（5,2）波形得到（5,3）波形，經過IIR濾波器后得到（5,4）時域波形。

圖6 matlab頻域分析如圖6所示，（6,3）與（6,4）相比3khz大概被削弱25DB左右。

圖7 IIR幅頻響應分析

如圖7所示，通過幅頻響應（7,2）可知3khz所在位置大概被削弱25DB。

4 FPGA設計驗證

FPGA IIR.v設計源碼：

`timescale 1ps/1ps
module iir(
input mclk,//45.1584MHZ
input reset_n,
input signed[31:0] pcm_in,
output signed[31:0] pcm_out
);

localparam LAST_CYCLE = 1023;
reg [9:0] i;

wire signed [15:0] b1,b2,b3,b4;
wire signed [15:0] a2,a3,a4;

wire signed [31:0] xn;
reg signed [31:0] xn1,xn2,xn3;
reg signed [31:0] yn,yn1,yn2,yn3;
reg signed [47:0] r_x1;
reg signed [47:0] r_x2;
reg signed [47:0] r_y;
reg signed [47:0] r_s;
reg signed [47:0] r_s1;

//coffe b
assign b1 = 3;
assign b2 = 8;
assign b3 = 8;
assign b4 = 3;

//coffe a
assign a2 = -22243;
assign a3 = 20231;
assign a4 = -6159;

assign xn = pcm_in;
assign pcm_out = yn;


always @(posedge mclk or negedge reset_n) begin
if(reset_n == 1'b0) begin
i <= 0;
xn1 <= 0;
xn2 <= 0;
xn3 <= 0;

yn <= 0;
yn1 <= 0;
yn2 <= 0;
yn3 <= 0;

r_x1 <= 0;
r_x2 <= 0;
r_y <= 0;
r_s <= 0;
r_s1 <= 0;
end
else begin
i<= i 1;
if(i==1) begin
r_x1 <= b1*(xn xn3);
r_x2 <= b2*(xn1 xn2);//Zero(n)
r_y <= a2*yn1 a3*yn2 a4*yn3;//Pole(n)
$display("r_x1 = %d,r_x2 = %d,r_y = %d",r_x1,r_x2,r_y);
end
if(i==2) begin
r_s <= r_x1 r_x2 - r_y;//8192y(n) = Zero(n)-Pole(n)
$display("%d",r_s);
end
if(i==3) r_s1 <= (r_s>>13);
if(i==4) yn <= r_s1[31:0];
if(i==5) begin //pipeline
xn1 <= xn;
xn2 <= xn1;
xn3 <= xn2;

yn1 <= yn;
yn2 <= yn1;
yn3 <= yn2;
end
end
end
Endmodule

實驗結果：
Modelsim波

圖8 modelsim 時域波形由圖8可知，pcm_out1(1khz 3khz)經過IIR濾波器后濾除了3khz，設計成功。
Matlab結果分析：

圖9 FPGA輸入波形matlab時域分析

圖10 FPGA結果數(shù)據(jù)分析由圖10的頻域分析結果可知3khz大概被削弱25db左右，設計成功。

感謝閱讀，別走！點贊、關注、轉發(fā)后再走吧