失眠网 > FIR内插滤波器的FPGA实现(一)-matlab实现

FIR内插滤波器的FPGA实现(一)-matlab实现

时间：2022-11-11 14:37:10

一、FIR内插滤波器

FIR内插滤波器是一种基本的插值方法，主要有两个步骤：

1）在输入的每一个初始采样中间插入L个零点；

2）插零后的数据经过低通滤波器。

二、matlab代码实现

1）插值前的信号：

T=0.5; %信号持续时间Fs=2000; %采样频率N=T*Fs; %采样点数x=(0:N-1)/Fs;f=200; %信号频率y=sin(2*pi*f*x);f=(0:T*Fs-1)/T-Fs/2;fft_y=2*abs(fftshift(fft(y)))/N;figure(1)subplot(211)plot(x,y,x,y,'c*');xlim([0,0.05]);subplot(212)plot(f,fft_y);

运行结果如图(时域、频域)：

我们的目标是把采样频率提升五十倍。

2）插零点：

%% 插0点L=49;%插值的点数y_1=[];for i=1:length(y)y_1=[y_1,y(i),zeros(1,L)];endfft_y1=2*abs(fftshift(fft(y_1)))/N/(L+1);f_1=(0:T*Fs*(L+1)-1)/T-Fs*(L+1)/2;figure(2)plot(f_1,fft_y1);

得到的频域图：

可以看出内插零点在频域会实现频谱扩展，因此也可以解释为什么低通滤波器可以还原出信号。

3）滤波

使用matlab设计滤波器，在命令框输入filterdesigner即可。

采样频率为100khz，通带、截止频率为1khz和2khz，滤波器形式为直接型，设计出的滤波器为254阶次的。

y_2 = (L+1)*filter(fir,1,y_1);x_1=(0:N*(L+1)-1)/Fs/(L+1);fft_y3=2*abs(fftshift(fft(y_2)))/N/(L+1);figure(3)subplot(211)plot(x_1,y_2);xlim([0,0.05]);subplot(212)plot(f_1,fft_y3);

完成了低通滤波，信号的形状不变，但是采样频率增加了50倍，实现了升采样。

原信号(y)，低通滤波后的信号(y_2)如图，采样频率提升50倍：

三、基于硬件实现的考虑

1）从FIR滤波器公式原理分析

FIR滤波器的过程实际上是乘系数再相加的过程，M为滤波器的阶数，b(k)为滤波器系数，如图：

y ( n ) = ∑ k = 0 M b ( k ) x ( n − k ) y(n) =\sum_{k=0}^{M} b(k)x(n-k) y(n)=k=0∑Mb(k)x(n−k)

基于以上公式，我们可以写FIR滤波器的代码：

%% 卷积计算for ii=1:254y_4(ii)=(L+1)*sum(fir(1:ii).*fliplr(y_1(1:ii)));endfor ii=255:(L+1)*Ny_4(ii)=(L+1)*(sum(fir.*fliplr(y_1(ii-253:ii))));endfft_y4=2*abs(fftshift(fft(y_4)))/N/(L+1);figure(4)subplot(211)plot(x_1,y_4);xlim([0,0.05]);subplot(212)plot(f_1,fft_y4);

运行结果是完全一样的。注意，曲线一开始的变形(红框处，即使用matlab的fir滤波器函数也有)，是由于一开始输入的信号个数比滤波器的阶次小，多余的滤波器系数乘0导致的。

但是这样的滤波器消耗了极大的乘法器资源，因为阶次有254阶，即使在FPGA用FIR滤波器的IP，消耗的资源也是我们无法接受的。

2）简化的方法

我们注意到：由于内插是零点，因此绝大多数的乘法是不必要的，实际上对于这样一个滤波器，每次的有效乘法次数只有5-6次(254阶/插值倍数(50))；

因此我们可以做如下的改进：

1.将fir滤波器的阶次扩展为300阶，通过补零实现，这样这个滤波器每次需要使用六次乘法。

2.通过六次乘法和相加实现完全相等的结果。

代码：

fir_1=[fir,zeros(1,300-254)];for ii=301:(L+1)*Na=ceil((ii-1)/50);y_3(ii)=(L+1)*sum(fir_1(ii-(a-1)*50-1)*y(a)+fir_1(ii-(a-2)*50-1)*y(a-1)+fir_1(ii-(a-3)*50-1)*y(a-2)+fir_1(ii-(a-4)*50-1)*y(a-3)+fir_1(ii-(a-5)*50-1)*y(a-4)+fir_1(ii-(a-6)*50-1)*y(a-5));endfigure(5)plot(x_1,y_3);xlim([0,0.05]);

具体结果如下:

实现的结果除了前300个点，与matlab的滤波器完全一致，前300个点不一样的原因在本文的三(1) 已经说明了，具体的改进方法也可以参考三(1) 中的代码，这里没有给出。实际应用中，这个偏差也是影响很小的，因为在有用的信号之前往往也会有一些值作为滤波器的输入。