matlab：FIR（有限脉冲响应）滤波器

 1、FIR滤波器的简单介绍

本文只讲代码，不详细陈述原理。

简单的说，就是：将一个包含噪声的时域信号，在其频域上去除噪声，然后再回到时域，得到想要的原始信号，通过fft(傅里叶变化函数)、ifft(傅里叶逆变换函数)实现。

再具体一点就是：将包含噪声的时域信号进行傅里叶变换，得到频谱图，然后自己设计一个滤波器，去除频谱图中不需要的频率点，只留下原始信号的频率点，然后再进行傅里叶逆变换，返回时域。

接下来通过一个例子进行理解

2、示例

①假设现在有一个频率为1Hz、幅值为10的时域信号，叠加上一个高斯白噪声，代码如下

fs = 50; % 采样频率
t = 0:1/fs:1-1/fs;
s = 10*sin(2*pi*t); % 原信号
white_noise = wgn(length(t), 1, 0)'; % 白噪声 
y = s + white_noise; % 原信号＋白噪声

figure(1);
subplot(3,1,1)
plot(t,s);
xlabel('t/s');
title('原信号')
subplot(3,1,2)
plot(t,white_noise);
xlabel('t/s');
title('高斯白噪声')
subplot(3,1,3)
plot(t,y);
xlabel('t/s');
title('原信号叠加高斯白噪声')

结果如下图所示

②对该进行傅里叶变换，然后设计一个滤波器，滤除其他频点（即将其他频率的值归零），只留下频率为1Hz处的值，代码如下

N = 50; % 采样点数
f_x = (0:(N-1))*fs/N; % 频谱横轴坐标，其中fs/N为频率分辨率
S_f = fft(y(1:N));
figure(2);
subplot(3,1,1)
stem(f_x, abs(S_f), 'filled');
xlabel('F/Hz');
title('频谱（原信号叠加高斯白噪声）')
S_f_FIR = zeros(size(S_f)); % FIR滤波器的频谱
S_f_FIR(2) = 1;  % 留下频率为1Hz处的值
S_f_FIR(50) = 1;  % 频谱上与1Hz对称的点
subplot(3,1,2)
stem(f_x, S_f_FIR, 'filled');
xlabel('F/Hz');
title('频谱（FIR滤波器）')
subplot(3,1,3)
stem(f_x, abs(S_f).*S_f_FIR, 'filled');
xlabel('F/Hz');
title('频谱（信号经过FIR滤波）')

结果如下图所示

需要注意的是：在之前的“matlab傅里叶变换”的文章中提到过，频率分辨率（fs/N）需要小于等于所求频率，才能避免出现频率混叠的现象，故此处fs/N应该小于等于1Hz

③进行傅里叶逆变换，回到时域，代码如下

s_FIR = ifft(S_f_FIR,length(t+1))
figure(3)
subplot(3,1,1)
plot(t,s);
xlabel('t/s');
title('原信号')
subplot(3,1,2)
plot(t,s_FIR)
xlabel('t/s');
title('FIR滤波器的时域信号')
s_t_FIR = ifft(S_f.*S_f_FIR,length(t+1)) % 该处不能用abs(S_f)，否则会丢失相位，导致时域信号时移
subplot(3,1,3)
plot(t,s_t_FIR)
xlabel('t/s');
title('通过FIR滤波器所恢复的时域信号')

结果如下图所示

需要注意的是：进行傅里叶逆变换的时候，对于原信号的频谱，不能取其模（即abs），因为取模的话，会丢失时域信号的相位，导致时域信号发生时移